Aplicando regresión univariante, modelo saturado, regresión con selección de modelos y KNN y árbol en modo regresión a house_prices.csv
El objetivo de esta práctica es predecir la variable continua SalePrice a través de distintos métodos de regresión y de algoritmos en modo regresión.
Necesitaremos los siguientes paquetes:
# Borramos
rm(list = ls())
# Paquetes
library(skimr) # Resumen numérico
library(tidymodels) # Depuración datos
library(tidyverse) # Modelos
library(outliers) # Outliers
library(themis) # Sobremuestreo
library(parallel) # Fase de validación
library(doParallel) # Fase de validación
library(rpart.plot) # Visualización de árboles
library(performance)
library(ggthemes)
library(glue)
library(vip)
library(ggrepel)
library(Rmisc)
Los datos que usaremos provienen de los dataset house_prices_train.csv y house_prices_test.csv. Para el análisis exploratorio, bindearemos ambos dataset para poder analizar el total de las observaciones.
# Cargamos ambas particiones
house_train <- read_csv(file = "/Users/leztin/Desktop/DATOS/house_prices_train.csv")
house_test <- read_csv(file = "/Users/leztin/Desktop/DATOS/house_prices_test.csv")
# Creamos en test una columna con la variable objetivo `SalePrice`
house_test$SalePrice <- NA
# Bindeamos las particiones en un archivo completo
house_complete <- rbind(house_train, house_test)
Antes de tomar cualquier decisión con los datos, lo primero que haremos será echar un vistazo numérico a cómo se comportan las variables. Dado que dos de los métodos que vamos a utilizar para predecir son el árbol y el KNN en modo regresión, comprobaremos cómo se relacionan nuestras predictoras cuantitativas y cualitativas de cara a la creación de la receta y con tal de poder recategorizarlas.
Lo primero es conocer nuestras variables.
glimpse(house_complete)
Rows: 2,919
Columns: 52
$ Id <dbl> 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14,…
$ MSSubClass <dbl> 60, 20, 60, 70, 60, 50, 20, 60, 50, 190, 20, 6…
$ MSZoning <chr> "RL", "RL", "RL", "RL", "RL", "RL", "RL", "RL"…
$ LotFrontage <dbl> 65, 80, 68, 60, 84, 85, 75, NA, 51, 50, 70, 85…
$ LotArea <dbl> 8450, 9600, 11250, 9550, 14260, 14115, 10084, …
$ Street <chr> "Pave", "Pave", "Pave", "Pave", "Pave", "Pave"…
$ Alley <chr> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA…
$ LandContour <chr> "Lvl", "Lvl", "Lvl", "Lvl", "Lvl", "Lvl", "Lvl…
$ Utilities <chr> "AllPub", "AllPub", "AllPub", "AllPub", "AllPu…
$ LandSlope <chr> "Gtl", "Gtl", "Gtl", "Gtl", "Gtl", "Gtl", "Gtl…
$ Neighborhood <chr> "CollgCr", "Veenker", "CollgCr", "Crawfor", "N…
$ BldgType <chr> "1Fam", "1Fam", "1Fam", "1Fam", "1Fam", "1Fam"…
$ HouseStyle <chr> "2Story", "1Story", "2Story", "2Story", "2Stor…
$ OverallCond <dbl> 5, 8, 5, 5, 5, 5, 5, 6, 5, 6, 5, 5, 6, 5, 5, 8…
$ YearBuilt <dbl> 2003, 1976, 2001, 1915, 2000, 1993, 2004, 1973…
$ YearRemodAdd <dbl> 2003, 1976, 2002, 1970, 2000, 1995, 2005, 1973…
$ RoofStyle <chr> "Gable", "Gable", "Gable", "Gable", "Gable", "…
$ RoofMatl <chr> "CompShg", "CompShg", "CompShg", "CompShg", "C…
$ Exterior1st <chr> "VinylSd", "MetalSd", "VinylSd", "Wd Sdng", "V…
$ ExterCond <chr> "TA", "TA", "TA", "TA", "TA", "TA", "TA", "TA"…
$ BsmtQual <chr> "Gd", "Gd", "Gd", "TA", "Gd", "Gd", "Ex", "Gd"…
$ BsmtCond <chr> "TA", "TA", "TA", "Gd", "TA", "TA", "TA", "TA"…
$ TotalBsmtSF <dbl> 856, 1262, 920, 756, 1145, 796, 1686, 1107, 95…
$ Heating <chr> "GasA", "GasA", "GasA", "GasA", "GasA", "GasA"…
$ HeatingQC <chr> "Ex", "Ex", "Ex", "Gd", "Ex", "Ex", "Ex", "Ex"…
$ CentralAir <chr> "Y", "Y", "Y", "Y", "Y", "Y", "Y", "Y", "Y", "…
$ Electrical <chr> "SBrkr", "SBrkr", "SBrkr", "SBrkr", "SBrkr", "…
$ `1stFlrSF` <dbl> 856, 1262, 920, 961, 1145, 796, 1694, 1107, 10…
$ `2ndFlrSF` <dbl> 854, 0, 866, 756, 1053, 566, 0, 983, 752, 0, 0…
$ GrLivArea <dbl> 1710, 1262, 1786, 1717, 2198, 1362, 1694, 2090…
$ FullBath <dbl> 2, 2, 2, 1, 2, 1, 2, 2, 2, 1, 1, 3, 1, 2, 1, 1…
$ HalfBath <dbl> 1, 0, 1, 0, 1, 1, 0, 1, 0, 0, 0, 0, 0, 0, 1, 0…
$ BedroomAbvGr <dbl> 3, 3, 3, 3, 4, 1, 3, 3, 2, 2, 3, 4, 2, 3, 2, 2…
$ KitchenQual <chr> "Gd", "TA", "Gd", "Gd", "Gd", "TA", "Gd", "TA"…
$ TotRmsAbvGrd <dbl> 8, 6, 6, 7, 9, 5, 7, 7, 8, 5, 5, 11, 4, 7, 5, …
$ Functional <chr> "Typ", "Typ", "Typ", "Typ", "Typ", "Typ", "Typ…
$ Fireplaces <dbl> 0, 1, 1, 1, 1, 0, 1, 2, 2, 2, 0, 2, 0, 1, 1, 0…
$ FireplaceQu <chr> NA, "TA", "TA", "Gd", "TA", NA, "Gd", "TA", "T…
$ GarageYrBlt <dbl> 2003, 1976, 2001, 1998, 2000, 1993, 2004, 1973…
$ GarageCars <dbl> 2, 2, 2, 3, 3, 2, 2, 2, 2, 1, 1, 3, 1, 3, 1, 2…
$ GarageArea <dbl> 548, 460, 608, 642, 836, 480, 636, 484, 468, 2…
$ GarageQual <chr> "TA", "TA", "TA", "TA", "TA", "TA", "TA", "TA"…
$ PavedDrive <chr> "Y", "Y", "Y", "Y", "Y", "Y", "Y", "Y", "Y", "…
$ OpenPorchSF <dbl> 61, 0, 42, 35, 84, 30, 57, 204, 0, 4, 0, 21, 0…
$ PoolArea <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
$ PoolQC <chr> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA…
$ Fence <chr> NA, NA, NA, NA, NA, "MnPrv", NA, NA, NA, NA, N…
$ MoSold <dbl> 2, 5, 9, 2, 12, 10, 8, 11, 4, 1, 2, 7, 9, 8, 5…
$ YrSold <dbl> 2008, 2007, 2008, 2006, 2008, 2009, 2007, 2009…
$ SaleType <chr> "WD", "WD", "WD", "WD", "WD", "WD", "WD", "WD"…
$ SaleCondition <chr> "Normal", "Normal", "Normal", "Abnorml", "Norm…
$ SalePrice <dbl> 208500, 181500, 223500, 140000, 250000, 143000…
| Variable | Significado | Variable | Significado |
|---|---|---|---|
SalePrice |
precio de venta de la propiedad en dólares (variable objetivo) | Heating |
tipo de calefacción |
MSSubClass |
clase de construcción | HeatingQC |
calidad y condición de la calefacción |
MSZoning |
clasificación general de zonificación | CentralAir |
aire acondicionado central |
LotFrontage |
pies lineales de calle conectados a la propiedad | Electrical |
sistema eléctrico |
LotArea |
tamaño de la propiedad en pies cuadrados | 1stFlrSF |
pies cuadrados del primer piso |
Street |
tipo de acceso a la carretera | 2ndFlrSF |
pies cuadrados del segundo piso |
Alley |
tipo de acceso al callejón | LowQualFinSF |
pies cuadrados terminados de baja calidad (todos los pisos) |
LotShape |
forma general de la propiedad | GrLivArea |
pies cuadrados de superficie habitable sobre el nivel del suelo (suelo) |
LandContour |
planitud de la propiedad | BsmtFullBath |
baños completos del sótano |
Utilities |
tipos de utilidades disponibles | BsmtHalfBath |
medios baños del sótano |
LotConfig |
configuración de las propiedades | FullBath |
baños completos sobre el nivel del suelo |
LandSlope |
pendiente de la propiedad | HalfBath |
medios baños sobre el nivel del suelo |
Neighborhood |
ubicaciones físicas dentro de los límites de la ciudad de Ames | Bedroom |
número de dormitorios por encima del nivel del sótano |
Condition1 |
proximidad a carretera principal o vía férrea | GarageFinish |
acabado interior del garaje |
Condition2 |
proximidad a la carretera principal o vía férrea (si hay una segunda presente) | GarageCars |
tamaño del garaje en capacidad de automóviles |
BldgType |
tipo de vivienda | GarageArea |
tamaño del garaje en pies cuadrados |
HouseStyle |
estilo de vivienda | GarageQual |
calidad de garaje |
OverallQual |
calidad general del material y del acabado | GarageCond |
estado del garaje |
OverallCond |
calificación de las condiciones generales | PavedDrive |
entrada pavimentada |
YearBuilt |
fecha de construcción original | WoodDeckSF |
área de la cubierta de madera en pies cuadrados |
YearRemodAdd |
fecha de remodelación | OpenPorchSF |
área de porche abierto en pies cuadrados |
RoofStyle |
tipo de techo | EnclosedPorch |
área de porche cerrado en pies cuadrados |
RoofMatl |
material del techo | 3SsnPorch |
área de porche de tres estaciones en pies cuadrados |
Exterior1st |
revestimiento exterior de la casa | ScreenPorch |
área de porche de pantalla en pies cuadrados |
Exterior2nd |
revestimiento exterior de la casa (si hay más de un material) | PoolArea |
área de la piscina en pies cuadrados |
MasVnrType |
tipo de chapa de mampostería | PoolQC |
calidad de la piscina |
MasVnrArea |
área de revestimiento de mampostería en pies cuadrados | Fence |
calidad de la cerca |
ExterQual |
calidad del material exterior | MiscFeature |
miscelánea no incluida en otras categorías |
ExterCond |
estado actual del material en el exterior | MiscVal |
miscelánea referente al valor |
Foundation |
tipo de cimentación | MoSold |
mes de venta |
BsmtQual |
altura del sótano | YrSold |
año de venta |
BsmtCond |
estado general del sótano | SaleType |
tipo de venta |
BsmtExposure |
paredes de los sótanos a nivel del jardín o de la salida | SaleCondition |
condiciones de venta |
BsmtFinType1 |
calidad del área terminada del sótano | Kitchen |
número de cocinas |
BsmtFinSF1 |
tipo 1: pies cuadrados terminados | KitchenQual |
calidad de la cocina |
BsmtFinType2 |
calidad de la segunda área terminada (si está presente) | TotRmsAbvGrd |
total de habitaciones sobre rasante (no incluye baños) |
BsmtFinSF2 |
tipo 2: pies cuadrados terminados | Functional |
calificación de la funcionalidad del hogar |
BsmtUnfSF |
pies cuadrados sin terminar del área del sótano | Fireplaces |
número de chimeneas |
TotalBsmtSF |
pies cuadrados totales del área del sótano | FireplaceQu |
calidad de chimenea |
GarageType |
ubicación del garaje | GarageYrBlt |
año en que se construyó el garaje |
El objetivo será predecir el precio de una vivienda en función de sus características, por lo que SalesPrice será nuestra variable objetivo.
En primer lugar comprobaremos cómo se distribuyen los valores de la objetivo.
# Objetivo: predecir el precio de una vivienda en función de sus características
house_complete |>
ggplot(aes(x = SalePrice)) +
geom_density(alpha = .8, fill="#EB9891") +
labs(title = "Distribución del precio de las viviendas", x = "Precio", y = NULL) +
theme_minimal() +
theme(text = element_text(face = "bold", size = 15), plot.title = element_text(hjust = 0.5),
axis.title.x = element_text(vjust=-0.5)) +
scale_x_continuous(labels = comma_format(big.mark = " ", decimal.mark = ",")) +
scale_y_continuous(labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_vline(aes(xintercept = mean(SalePrice, na.rm = T), linetype = "Media"), colour = "black", size = .8) +
geom_vline(aes(xintercept = median(SalePrice, na.rm = T), linetype = "Mediana"), colour = "black", size = .8) +
scale_linetype_manual(name = "Medidas", values = c(Media = "solid", Mediana = "dotted"))
En este caso, como la distribución parece bastante asimétrica, la mediana será la medida más representativa. Esta se encuentra en torno a los 160 000 dólares por vivienda.
Tras esta pequeña aproximación al dataset, comienza la primera fase de la metodología SEMMA para el depurado de nuestros datos. En este primer apartado observaremos a grosso modo si existen problemas de codificación en el dataset. Lo que más llama la atención con respecto a la codificación de las modalidades de ciertas variables es el gran número de ausentes que presenta el dataset.
ausentes <-
apply(house_complete, 2, function(x) sum(is.na(x)))
ausentes_tb <-
tibble(Variable = names(house_complete), Ausentes = ausentes) |>
filter(Ausentes > 0)
ausentes_tb
# A tibble: 20 × 2
Variable Ausentes
<chr> <int>
1 MSZoning 4
2 LotFrontage 486
3 Alley 2721
4 Utilities 2
5 Exterior1st 1
6 BsmtQual 81
7 BsmtCond 82
8 TotalBsmtSF 1
9 Electrical 1
10 KitchenQual 1
11 Functional 2
12 FireplaceQu 1420
13 GarageYrBlt 159
14 GarageCars 1
15 GarageArea 1
16 GarageQual 159
17 PoolQC 2909
18 Fence 2348
19 SaleType 1
20 SalePrice 1459
ausentes_tb |>
filter(Ausentes > 10) |>
ggplot(aes(x = Variable, y = Ausentes)) +
geom_col(position = "dodge", fill="#EB9891", color = "black") +
labs(title = "> 10 valores ausentes por variable", x = "Variable",
y = "Cantidad de valores ausentes") +
theme_minimal() +
theme(text = element_text(face = "bold", size = 9), plot.title = element_text(hjust = 0.5),
axis.title.x = element_text(vjust=-0.5), axis.title.y = element_text(vjust=2)) +
scale_x_discrete(guide = guide_axis(n.dodge = 2))
Como se puede observar, 20 del total de variables presenta algún dato ausente. La mitad de esos 20 presenta más de 10 ausentes entre sus categorías. Ante esta cantidad de valores vacíos, debemos consultar primero si se tratan de categorías reales mal codificadas o ausentes reales.
Según el documento explicativo, las variables Alley, BsmtQual, BsmtCond, FireplaceQu, GarageQual, PoolQC y Fence tienen codificadas como valores ausentes categorías que realmente no lo son. Vamos a modificar esta cuestión.
house <-
house_complete |>
mutate(Alley =
ifelse(is.na(Alley), "None", Alley)) |>
mutate(BsmtQual =
ifelse(is.na(BsmtQual), "None", BsmtQual)) |>
mutate(BsmtCond =
ifelse(is.na(BsmtCond), "None", BsmtCond)) |>
mutate(FireplaceQu =
ifelse(is.na(FireplaceQu), "None", FireplaceQu)) |>
mutate(GarageQual =
ifelse(is.na(GarageQual), "None", GarageQual)) |>
mutate(PoolQC =
ifelse(is.na(PoolQC), "None", PoolQC)) |>
mutate(Fence =
ifelse(is.na(Fence), "None", Fence))
ausentes <-
apply(house, 2, function(x) sum(is.na(x)))
ausentes_tb <-
tibble(Variable = names(house), Ausentes = ausentes) |>
filter(Ausentes > 0)
ausentes_tb
# A tibble: 13 × 2
Variable Ausentes
<chr> <int>
1 MSZoning 4
2 LotFrontage 486
3 Utilities 2
4 Exterior1st 1
5 TotalBsmtSF 1
6 Electrical 1
7 KitchenQual 1
8 Functional 2
9 GarageYrBlt 159
10 GarageCars 1
11 GarageArea 1
12 SaleType 1
13 SalePrice 1459
Ahora sí se ha reducido el número de variables con ausentes de 20 a 13. Al resto, le imputaremos lo que le corresponda (media, moda o mediana) en la fase pertinente.
Tras ello, comprobaremos que todas las variables estén codificadas en su tipología correcta: debemos decidir si las variables de tipo texto son realmente variables cualitativas (factores).
Rows: 2,919
Columns: 28
$ MSZoning <chr> "RL", "RL", "RL", "RL", "RL", "RL", "RL", "RL"…
$ Street <chr> "Pave", "Pave", "Pave", "Pave", "Pave", "Pave"…
$ Alley <chr> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA…
$ LandContour <chr> "Lvl", "Lvl", "Lvl", "Lvl", "Lvl", "Lvl", "Lvl…
$ Utilities <chr> "AllPub", "AllPub", "AllPub", "AllPub", "AllPu…
$ LandSlope <chr> "Gtl", "Gtl", "Gtl", "Gtl", "Gtl", "Gtl", "Gtl…
$ Neighborhood <chr> "CollgCr", "Veenker", "CollgCr", "Crawfor", "N…
$ BldgType <chr> "1Fam", "1Fam", "1Fam", "1Fam", "1Fam", "1Fam"…
$ HouseStyle <chr> "2Story", "1Story", "2Story", "2Story", "2Stor…
$ RoofStyle <chr> "Gable", "Gable", "Gable", "Gable", "Gable", "…
$ RoofMatl <chr> "CompShg", "CompShg", "CompShg", "CompShg", "C…
$ Exterior1st <chr> "VinylSd", "MetalSd", "VinylSd", "Wd Sdng", "V…
$ ExterCond <chr> "TA", "TA", "TA", "TA", "TA", "TA", "TA", "TA"…
$ BsmtQual <chr> "Gd", "Gd", "Gd", "TA", "Gd", "Gd", "Ex", "Gd"…
$ BsmtCond <chr> "TA", "TA", "TA", "Gd", "TA", "TA", "TA", "TA"…
$ Heating <chr> "GasA", "GasA", "GasA", "GasA", "GasA", "GasA"…
$ HeatingQC <chr> "Ex", "Ex", "Ex", "Gd", "Ex", "Ex", "Ex", "Ex"…
$ CentralAir <chr> "Y", "Y", "Y", "Y", "Y", "Y", "Y", "Y", "Y", "…
$ Electrical <chr> "SBrkr", "SBrkr", "SBrkr", "SBrkr", "SBrkr", "…
$ KitchenQual <chr> "Gd", "TA", "Gd", "Gd", "Gd", "TA", "Gd", "TA"…
$ Functional <chr> "Typ", "Typ", "Typ", "Typ", "Typ", "Typ", "Typ…
$ FireplaceQu <chr> NA, "TA", "TA", "Gd", "TA", NA, "Gd", "TA", "T…
$ GarageQual <chr> "TA", "TA", "TA", "TA", "TA", "TA", "TA", "TA"…
$ PavedDrive <chr> "Y", "Y", "Y", "Y", "Y", "Y", "Y", "Y", "Y", "…
$ PoolQC <chr> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA…
$ Fence <chr> NA, NA, NA, NA, NA, "MnPrv", NA, NA, NA, NA, N…
$ SaleType <chr> "WD", "WD", "WD", "WD", "WD", "WD", "WD", "WD"…
$ SaleCondition <chr> "Normal", "Normal", "Normal", "Abnorml", "Norm…
En el dataset encontramos varias variables divididas en modalidades del tipo Ex (excelente) > Gd (bueno) > TA (normal) > Fa (regular) > Po (pésimo). En concreto, se tratan de las variables ExterCond, BsmtQual, BsmtCond, HeatingQC, KitchenQual, FireplaceQu, GarageQual y PoolQC. Para ampliar el análisis de correlaciones posterior para con nuestra variable objetivo, lo que haremos será transformar todas estas variables cualitativas en ordinales. Como la mayoría sigue la misma estructura, crearemos un vector que asigne a cada modalidad un número del 0 al 5, siendo 0 igual a nada, y siendo 5 igual a excelente:
quality <-
c("None" = 0, "Po" = 1, "Fa" = 2, "TA" = 3, "Gd" = 4, "Ex" = 5)
Aplicamos el vector con revalue() a cada una de las variables:
house$ExterCond <-
as.integer(revalue(house$ExterCond, quality))
house$BsmtQual <-
as.integer(revalue(house$BsmtQual, quality))
house$BsmtCond <-
as.integer(revalue(house$BsmtCond, quality))
house$HeatingQC <-
as.integer(revalue(house$HeatingQC, quality))
house$KitchenQual <-
as.integer(revalue(house$KitchenQual, quality))
house$FireplaceQu <-
as.integer(revalue(house$FireplaceQu, quality))
house$GarageQual <-
as.integer(revalue(house$GarageQual, quality))
house$PoolQC <-
as.integer(revalue(house$PoolQC, quality))
Para el resto de variables no ordinales las convertimos a factor.
house <-
house |>
mutate_if(~!is.numeric(.), as.factor)
Daremos un repaso también a las numéricas.
Rows: 2,919
Columns: 24
$ Id <dbl> 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, …
$ MSSubClass <dbl> 60, 20, 60, 70, 60, 50, 20, 60, 50, 190, 20, 60…
$ LotFrontage <dbl> 65, 80, 68, 60, 84, 85, 75, NA, 51, 50, 70, 85,…
$ LotArea <dbl> 8450, 9600, 11250, 9550, 14260, 14115, 10084, 1…
$ OverallCond <dbl> 5, 8, 5, 5, 5, 5, 5, 6, 5, 6, 5, 5, 6, 5, 5, 8,…
$ YearBuilt <dbl> 2003, 1976, 2001, 1915, 2000, 1993, 2004, 1973,…
$ YearRemodAdd <dbl> 2003, 1976, 2002, 1970, 2000, 1995, 2005, 1973,…
$ TotalBsmtSF <dbl> 856, 1262, 920, 756, 1145, 796, 1686, 1107, 952…
$ `1stFlrSF` <dbl> 856, 1262, 920, 961, 1145, 796, 1694, 1107, 102…
$ `2ndFlrSF` <dbl> 854, 0, 866, 756, 1053, 566, 0, 983, 752, 0, 0,…
$ GrLivArea <dbl> 1710, 1262, 1786, 1717, 2198, 1362, 1694, 2090,…
$ FullBath <dbl> 2, 2, 2, 1, 2, 1, 2, 2, 2, 1, 1, 3, 1, 2, 1, 1,…
$ HalfBath <dbl> 1, 0, 1, 0, 1, 1, 0, 1, 0, 0, 0, 0, 0, 0, 1, 0,…
$ BedroomAbvGr <dbl> 3, 3, 3, 3, 4, 1, 3, 3, 2, 2, 3, 4, 2, 3, 2, 2,…
$ TotRmsAbvGrd <dbl> 8, 6, 6, 7, 9, 5, 7, 7, 8, 5, 5, 11, 4, 7, 5, 5…
$ Fireplaces <dbl> 0, 1, 1, 1, 1, 0, 1, 2, 2, 2, 0, 2, 0, 1, 1, 0,…
$ GarageYrBlt <dbl> 2003, 1976, 2001, 1998, 2000, 1993, 2004, 1973,…
$ GarageCars <dbl> 2, 2, 2, 3, 3, 2, 2, 2, 2, 1, 1, 3, 1, 3, 1, 2,…
$ GarageArea <dbl> 548, 460, 608, 642, 836, 480, 636, 484, 468, 20…
$ OpenPorchSF <dbl> 61, 0, 42, 35, 84, 30, 57, 204, 0, 4, 0, 21, 0,…
$ PoolArea <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
$ MoSold <dbl> 2, 5, 9, 2, 12, 10, 8, 11, 4, 1, 2, 7, 9, 8, 5,…
$ YrSold <dbl> 2008, 2007, 2008, 2006, 2008, 2009, 2007, 2009,…
$ SalePrice <dbl> 208500, 181500, 223500, 140000, 250000, 143000,…
Nos llama la atención la variable MSSubClass, que identifica el tipo de vivienda objeto de venta. Aunque este codificada como numérica, realmente se trata de una variable cualitativa, por lo que vamos a convertirla.
class <-
c("20" = "1 story 1946+", "30" = "1 story 1945-",
"40" = "1 story fin attic", "45" = "1,5 story unf",
"50" = "1,5 story fin", "60" = "2 story 1946+",
"70" = "2 story 1945-", "75" = "2,5 story all ages",
"80" = "split/multi level", "85" = "split foyer",
"90" = "duplex all style/age", "120" = "1 story PUD 1946+",
"150" = "1,5 story PUD all", "160" = "2 story PUD 1946+",
"180" = "PUD multilevel", "190" = "2 family conversion")
house$MSSubClass <-
as.factor(house$MSSubClass)
house$MSSubClass <-
revalue(house$MSSubClass, class)
Ahora sí, todas las variables cualitativas y cuantitativas están bien codificadas.
Una vez asignado a cada variable su tipología correspondiente, pasaremos a analizar las variables cuantitativas del dataset. Se analizará ante todo cómo afecta cada variable a nuestra variable objetivo (SalePrice). Este análisis servirá, ante todo, para recategorizar las predictoras numéricas a la hora de crear la receta para el árbol en modo regresión.
Antes de nada, al tener que analizar un dataset con tantas variables, comprobaremos los posibles problemas de colinealidad entre las predictoras numéricas con tal de eliminar las que repitan información. También, al tener una variable continua como objetivo, comprobaremos cuáles de las numéricas tienen una mayor correlación con ella con tal de mantenerlas y analizarlas en profundidad.
library(corrplot)
cor_matrix |>
corrplot(method = "number", tl.cex = 0.55, number.cex = 0.7, type = "lower")
Como se puede observar, existen bastantes variables con correlaciones superiores a 0.49. Vamos a echarles un vistazo:
# Seleccionar solo las entradas de la matriz de correlación con valores mayores a 0.49
filtered_matrix <- abs(cor_matrix)
filtered_matrix[cor_matrix <= 0.49] <- NA
# Creamos el gráfico de correlación
corrplot(filtered_matrix, method = "number", type = "lower",tl.cex = 0.55, number.cex = 0.7, na.label=" ")
En primer lugar, a partir del gráfico podemos observar como las variables YearsBuilt, YearRemodAdd, BsmtQual, TotalBsmtSF, 1stFlrSF, GrLivArea, FullBath, KitchenQual, TotRmsAbvGrd, FireplaceQu GarageCars y GarageArea son las que mantienen correlaciones más fuertes (superiores a 0.5 en valor absoluto) con nuestra variable objetivo SalePrice. Podríamos añadir también las variables Fireplaces y GarageYrBlt ya que mantienen correlaciones de 0.47 y 0.49, respectivamente.
En segundo lugar, analizaremos las correlaciones entre las propias predictoras con tal de deshacernos de información redundante. Nos desharemos de aquellas variables que mantengan correlaciones superiores a 0.7 con otras variables. De un par de variables muy correlacionadas, eliminaremos la que menor correlación mantenga con nuestra objetivo (SalePrice).
YearsBuilt mantiene una correlación de 0.83 con GarageYrBlt. Ante este par, eliminaremos GarageYrBlt.TotalBsmtSF mantiene una correlación de 0.8 con 1stFlrSF. Ante este par, eliminaremos TotalBsmtSF.GrLivArea mantiene una correlación de 0.81 con TotRmsAbvGrd. Ante este par, eliminaremos TotRmsAbvGrd.Fireplaces mantiene una correlación de 0.86 con FireplaceQu. Ante este par, eliminaremos Fireplaces.GarageCars mantiene una correlación de 0.89 con GarageArea. Ante este par, eliminaremos GarageArea.PoolArea mantiene una correlación de 0.79 con PoolQC. Ante este par, eliminaremos PoolArea.house <-
house |>
select(-c(GarageYrBlt, TotalBsmtSF, TotRmsAbvGrd, Fireplaces, GarageArea, PoolArea))
house |>
select(where(is.numeric)) |>
glimpse()
Rows: 2,919
Columns: 25
$ Id <dbl> 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, …
$ LotFrontage <dbl> 65, 80, 68, 60, 84, 85, 75, NA, 51, 50, 70, 85,…
$ LotArea <dbl> 8450, 9600, 11250, 9550, 14260, 14115, 10084, 1…
$ OverallCond <dbl> 5, 8, 5, 5, 5, 5, 5, 6, 5, 6, 5, 5, 6, 5, 5, 8,…
$ YearBuilt <dbl> 2003, 1976, 2001, 1915, 2000, 1993, 2004, 1973,…
$ YearRemodAdd <dbl> 2003, 1976, 2002, 1970, 2000, 1995, 2005, 1973,…
$ ExterCond <int> 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3,…
$ BsmtQual <int> 4, 4, 4, 3, 4, 4, 5, 4, 3, 3, 3, 5, 3, 4, 3, 3,…
$ BsmtCond <int> 3, 3, 3, 4, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3,…
$ HeatingQC <int> 5, 5, 5, 4, 5, 5, 5, 5, 4, 5, 5, 5, 3, 5, 3, 5,…
$ `1stFlrSF` <dbl> 856, 1262, 920, 961, 1145, 796, 1694, 1107, 102…
$ `2ndFlrSF` <dbl> 854, 0, 866, 756, 1053, 566, 0, 983, 752, 0, 0,…
$ GrLivArea <dbl> 1710, 1262, 1786, 1717, 2198, 1362, 1694, 2090,…
$ FullBath <dbl> 2, 2, 2, 1, 2, 1, 2, 2, 2, 1, 1, 3, 1, 2, 1, 1,…
$ HalfBath <dbl> 1, 0, 1, 0, 1, 1, 0, 1, 0, 0, 0, 0, 0, 0, 1, 0,…
$ BedroomAbvGr <dbl> 3, 3, 3, 3, 4, 1, 3, 3, 2, 2, 3, 4, 2, 3, 2, 2,…
$ KitchenQual <int> 4, 3, 4, 4, 4, 3, 4, 3, 3, 3, 3, 5, 3, 4, 3, 3,…
$ FireplaceQu <int> 0, 3, 3, 4, 3, 0, 4, 3, 3, 3, 0, 4, 0, 4, 2, 0,…
$ GarageCars <dbl> 2, 2, 2, 3, 3, 2, 2, 2, 2, 1, 1, 3, 1, 3, 1, 2,…
$ GarageQual <int> 3, 3, 3, 3, 3, 3, 3, 3, 2, 4, 3, 3, 3, 3, 3, 3,…
$ OpenPorchSF <dbl> 61, 0, 42, 35, 84, 30, 57, 204, 0, 4, 0, 21, 0,…
$ PoolQC <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
$ MoSold <dbl> 2, 5, 9, 2, 12, 10, 8, 11, 4, 1, 2, 7, 9, 8, 5,…
$ YrSold <dbl> 2008, 2007, 2008, 2006, 2008, 2009, 2007, 2009,…
$ SalePrice <dbl> 208500, 181500, 223500, 140000, 250000, 143000,…
Como se puede observar, la variable Id es un id del número de viviendas que se incluyen en el dataset. Como no tiene interés alguno a fin de predecir nuestra variable objetivo, en la fase de modificación eliminaremos esta variable. La recuperaremos únicamente para cuando llegue el momento de subir las predicciones a Kaggle.
sum1 <- house |>
summarise(Variable = "GrLivArea", min_lead = min(GrLivArea), max_lead = max(GrLivArea))
sum2 <- house |>
summarise(Variable = "`1stFlrSF`", min_lead = min(`1stFlrSF`), max_lead = max(`1stFlrSF`))
sum3 <- house |>
summarise(Variable = "`2ndFlrSF`", min_lead = min(`2ndFlrSF`), max_lead = max(`2ndFlrSF`))
sum4 <- house |>
summarise(Variable = "LotArea", min_lead = min(LotArea), max_lead = max(LotArea))
sum5 <- house |>
drop_na(LotFrontage) |>
summarise(Variable = "LotFrontage", min_lead = min(LotFrontage), max_lead = max(LotFrontage))
sum6 <- house |>
summarise(Variable = "OpenPorchSF", min_lead = min (OpenPorchSF), max_lead = max(OpenPorchSF))
rbind(sum1, sum2, sum3, sum4, sum5, sum6)
Variable min_lead max_lead
1 GrLivArea 334 5642
2 `1stFlrSF` 334 5095
3 `2ndFlrSF` 0 2065
4 LotArea 1300 215245
5 LotFrontage 21 313
6 OpenPorchSF 0 742
Las seis variables miden áreas de ciertas zonas de la vivienda. Las variables 2ndFlrSF y OpenPorchSF toman el valor 0 cuando la vivienda no tiene segunda planta o porche. Graficaremos a continuación las seis variables para comprobar como se distribuyen en función de nuestra objetivo SalePrice.
a1 <- house[!is.na(house$SalePrice), ] |>
ggplot(aes(x = GrLivArea, y = SalePrice)) +
geom_point(col = "#EB9891") +
geom_smooth(method = "lm", se = FALSE, color = "black", aes(group = 1)) +
scale_y_continuous(breaks= seq(0, 800000, by = 100000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
scale_x_continuous(labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_text_repel(aes(label = ifelse(house$GrLivArea[!is.na(house$SalePrice)] > 4500,
rownames(house), ""))) +
theme_minimal()
a2 <- house[!is.na(house$SalePrice), ] |>
ggplot(aes(x = `1stFlrSF`, y = SalePrice)) +
geom_point(col = "#EB9891") +
geom_smooth(method = "lm", se = FALSE, color = "black", aes(group = 1)) +
scale_y_continuous(breaks= seq(0, 800000, by = 100000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
scale_x_continuous(labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_text_repel(aes(label = ifelse(house$`1stFlrSF`[!is.na(house$SalePrice)] > 4500,
rownames(house), ""))) +
theme_minimal()
a3 <- house[!is.na(house$SalePrice), ] |>
filter(`2ndFlrSF` > 0) |>
ggplot(aes(x = `2ndFlrSF`, y = SalePrice)) +
geom_point(col = "#EB9891") +
geom_smooth(method = "lm", se = FALSE, color = "black", aes(group = 1)) +
scale_y_continuous(breaks= seq(0, 800000, by = 100000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
scale_x_continuous(labels = comma_format(big.mark = " ", decimal.mark = ",")) +
theme_minimal()
a4 <- house[!is.na(house$SalePrice), ] |>
filter(LotArea < 100000) |>
ggplot(aes(x = LotArea, y = SalePrice)) +
geom_point(col = "#EB9891") +
geom_smooth(method = "lm", se = FALSE, color = "black", aes(group = 1)) +
scale_y_continuous(breaks= seq(0, 800000, by = 100000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
scale_x_continuous(labels = comma_format(big.mark = " ", decimal.mark = ",")) +
theme_minimal()
a5 <- house[!is.na(house$SalePrice), ] |>
ggplot(aes(x = LotFrontage, y = SalePrice)) +
geom_point(col = "#EB9891") +
geom_smooth(method = "lm", se = FALSE, color = "black", aes(group = 1)) +
scale_y_continuous(breaks= seq(0, 800000, by = 100000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
scale_x_continuous(labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_text_repel(aes(label = ifelse(house$LotFrontage[!is.na(house$SalePrice)] > 250,
rownames(house), ""))) +
theme_minimal()
a6 <- house[!is.na(house$SalePrice), ] |>
filter(OpenPorchSF > 0) |>
ggplot(aes(x = OpenPorchSF, y = SalePrice)) +
geom_point(col = "#EB9891") +
geom_smooth(method = "lm", se = FALSE, color = "black", aes(group = 1)) +
scale_y_continuous(breaks= seq(0, 800000, by = 100000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
scale_x_continuous(labels = comma_format(big.mark = " ", decimal.mark = ",")) +
theme_minimal()
multiplot(a1, a2, a3, a4, a5, a6, cols = 2)
Como se puede observar, no todas las variables evolucionan de la misma manera al relacionarlas con la variable SalePrice.
Para las variables GrLivArea, 1stFlrSF y 2ndFlrSF la correlación es muy directa: como es evidente, a más superficie habitable haya en la vivienda, mayor será su precio. De hecho, la variable GrLivArea parece ser resultado de la suma de las variables 1stFlrSF y 2ndFlrSF, que medían la superficie habitable de la vivienda en el segundo y el segundo piso (si hubiere). Para poder graficar correctamente la variable 2ndFlrSF, se han filtrado los valores iguales a 0 (lo que era equivalente a no tener segundo piso).
Para las variables LotArea y LotFrontage la relación también es directa, pero la evolución no es tan clara como en las anteriores variables. Es evidente que sigue siendo positiva porque ambas variables miden la superficie habitable y construida de la vivienda.
Por su parte, la variable OpenPorchSF no parece tener una influencia demasiado significativa sobre nuestra variable objetivo. A medida que el porche aumenta en pies cuadrados, el precio general de la vivienda no parece variar en exceso.
Como ya tratamos los ausentes de estas variables (a excepción de LotFrontage, a los que probablemente imputemos la mediana o la eliminemos por su elevado número de ausentes), dejaremos estas seis variables tal y como están, a excepción de OpenPorchSF, la cual seguramente eliminemos.
sum1 <- house |>
summarise(Variable = "OverallCond", min_lead = min(OverallCond), max_lead = max(OverallCond))
sum2 <- house |>
summarise(Variable = "ExterCond", min_lead = min(ExterCond), max_lead = max(ExterCond))
sum3 <- house |>
summarise(Variable = "BsmtQual", min_lead = min(BsmtQual), max_lead = max(BsmtQual))
sum4 <- house |>
summarise(Variable = "BsmtCond", min_lead = min(BsmtCond), max_lead = max(BsmtCond))
sum5 <- house |>
drop_na(KitchenQual) |>
summarise(Variable = "KitchenQual", min_lead = min(KitchenQual), max_lead = max(KitchenQual))
sum6 <- house |>
summarise(Variable = "GarageQual", min_lead = min(GarageQual), max_lead = max(GarageQual))
sum7 <- house |>
summarise(Variable = "FireplaceQu", min_lead = min(FireplaceQu), max_lead = max(FireplaceQu))
sum8 <- house |>
summarise(Variable = "PoolQC", min_lead = min(PoolQC), max_lead = max(PoolQC))
sum9 <- house |>
summarise(Variable = "HeatingQC", min_lead = min(HeatingQC), max_lead = max(HeatingQC))
rbind(sum1, sum2, sum3, sum4, sum5, sum6, sum7, sum8, sum9)
Variable min_lead max_lead
1 OverallCond 1 9
2 ExterCond 1 5
3 BsmtQual 0 5
4 BsmtCond 0 4
5 KitchenQual 2 5
6 GarageQual 0 5
7 FireplaceQu 0 5
8 PoolQC 0 5
9 HeatingQC 1 5
Las nueve variables miden las condiciones y calidades de ciertas zonas de la vivienda. Las variables BsmtQual, GarageQual, FireplaceQu y PoolQC toman el valor 0 cuando la vivienda no dispone de esas áreas. Graficaremos a continuación las nueve variables para comprobar como se distribuyen en función de nuestra objetivo SalePrice.
c1 <- house[!is.na(house$SalePrice), ] |>
ggplot(aes(x = factor(OverallCond), y = SalePrice)) +
geom_boxplot(col = "#EB9891") + labs(x = "Overall Condition") +
scale_y_continuous(breaks= seq(0, 800000, by = 100000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
theme_minimal()
c2 <- house[!is.na(house$SalePrice), ] |>
ggplot(aes(x = factor(ExterCond), y = SalePrice)) +
geom_boxplot(col = "#EB9891") + labs(x = "Exterior Condition") +
scale_y_continuous(breaks= seq(0, 800000, by = 100000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
theme_minimal()
c3 <- house[!is.na(house$SalePrice), ] |>
ggplot(aes(x = factor(BsmtQual), y = SalePrice)) +
geom_boxplot(col = "#EB9891") + labs(x = "Basement Quality") +
scale_y_continuous(breaks= seq(0, 800000, by = 100000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
theme_minimal()
c4 <- house[!is.na(house$SalePrice), ] |>
ggplot(aes(x = factor(BsmtCond), y = SalePrice)) +
geom_boxplot(col = "#EB9891") + labs(x = "Basement Condition") +
scale_y_continuous(breaks= seq(0, 800000, by = 100000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
theme_minimal()
c5 <- house[!is.na(house$SalePrice), ] |>
ggplot(aes(x = factor(KitchenQual), y = SalePrice)) +
geom_boxplot(col = "#EB9891") + labs(x = "Kitchen Quality") +
scale_y_continuous(breaks= seq(0, 800000, by = 100000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
theme_minimal()
c6 <- house[!is.na(house$SalePrice), ] |>
ggplot(aes(x = factor(GarageQual), y = SalePrice)) +
geom_boxplot(col = "#EB9891") + labs(x = "Garage Quality") +
scale_y_continuous(breaks= seq(0, 800000, by = 100000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
theme_minimal()
c7 <- house[!is.na(house$SalePrice), ] |>
ggplot(aes(x = factor(FireplaceQu), y = SalePrice)) +
geom_boxplot(col = "#EB9891") + labs(x = "Fireplace Quality") +
scale_y_continuous(breaks= seq(0, 800000, by=100000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
theme_minimal()
c8 <- house[!is.na(house$SalePrice), ] |>
ggplot(aes(x = factor(PoolQC), y = SalePrice)) +
geom_boxplot(col = "#EB9891") + labs(x = "Pool Quality") +
scale_y_continuous(breaks= seq(0, 800000, by=100000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
theme_minimal()
c9 <- house[!is.na(house$SalePrice), ] |>
ggplot(aes(x = factor(HeatingQC), y = SalePrice)) +
geom_boxplot(col = "#EB9891") + labs(x = "Heating Quality") +
scale_y_continuous(breaks= seq(0, 800000, by=100000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
theme_minimal()
multiplot(c1, c2, c3, c4, c5, c6, c7, c8, c9, cols = 3)
A priori, podríamos pensar que a mayor calidad en las zonas de la vivienda, mayor será el precio de la misma, pero no con todas parece cumplirse esta afirmación.
Para las variables BsmtQual, FireplaceQu y KitchenQual la correlación es muy directa y la progresión parece ligeramente ascendente: en estas tres variables se observa cierta correlación con la variable objetivo. De hecho, del grupo de variables que miden las calidades y las condiciones de la vivienda, estas tres eran las que mayor grado de correlación mantenían con SalePrice, de 0.59, 0.52 y 0.66, respectivamente.
Para las variables OverallCond, ExterCond y GarageQual la relación no parece tan evidente. En las tres parece que la relación comienza con una progresión ascendiente para las primeras categorías, pero llegados a un punto, esta se estanca. Para el caso de OverallCond, por ejemplo, mientras que las diferencias de precios respecto del rango 1 a 5 sí que parecen aumentar ligeramente, a partir del valor 5 la progresión se estanca e incluso disminuye para el resto de valores (6-10). Lo mismo ocurre para ExterCond con el valor 3, y a GarageQual con el 4. Además, estas modalidades presentan una gran cantidad de valores outliers.
La variable PoolQC por sí sola no parece ser de mucho interés. La mayoría de valores se concentran en 0, indicando con ello que en la mayoría de viviendas no hay piscina que evaluar. Se considerará posteriormente la posibilidad de sustituir las variables relativas a piscinas por una binaria que recoja únicamente la información de si se tiene o no con categorías Yes y No.
b1 <- ggplot(house, aes(x = factor(OverallCond))) +
geom_histogram(stat = "Count", fill = "#EB9891") +
labs(x="OverallCond") +
theme_minimal()
b2 <- ggplot(house, aes(x = factor(ExterCond))) +
geom_histogram(stat = "Count", fill = "#EB9891") +
labs(x="ExterCond") +
theme_minimal()
b3 <- ggplot(house, aes(x = factor(BsmtQual))) +
geom_histogram(stat = "Count", fill = "#EB9891") +
labs(x="BsmtQual") +
theme_minimal()
b4 <- ggplot(house, aes(x = factor(BsmtCond))) +
geom_histogram(stat = "Count", fill = "#EB9891") +
labs(x="BsmtCond") +
theme_minimal()
b5 <- ggplot(house, aes(x = factor(KitchenQual))) +
geom_histogram(stat = "Count", fill = "#EB9891") +
labs(x="KitchenQual") +
theme_minimal()
b6 <- ggplot(house, aes(x = factor(GarageQual))) +
geom_histogram(stat = "Count", fill = "#EB9891") +
labs(x="GarageQual") +
theme_minimal()
b7 <- ggplot(house, aes(x = factor(FireplaceQu))) +
geom_histogram(stat = "Count", fill = "#EB9891") +
labs(x="FireplaceQu") +
theme_minimal()
b8 <- ggplot(house, aes(x = factor(PoolQC))) +
geom_histogram(stat = "Count", fill = "#EB9891") +
labs(x="PoolQC") +
theme_minimal()
b9 <- ggplot(house, aes(x = factor(HeatingQC))) +
geom_histogram(stat = "Count", fill = "#EB9891") +
labs(x="HeatingQC") +
theme_minimal()
multiplot(b1, b2, b3, b4, b5, b6, b7, b8, b9, cols = 3)
Por otro lado, si se atiende a la representación de cada categoría sobre el total de registros, en las variables GarageQual, BsmtCond o ExterCond la mayor parte las gobiernan la modalidad 3. Posteriormente, en la fase de la preparación de la receta, se podrían agrupar algunos de estos niveles con tal de romper la concentración de los datos o, por el contrario, eliminar directamente este tipo de variables.
De momento, mantendremos todas las variables a excepción de PoolQC, que será sustituida por una binaria.
sum1 <- house |>
summarise(Variable = "YearRemodAdd", min_lead = min(YearRemodAdd), max_lead = max(YearRemodAdd))
sum2 <- house |>
summarise(Variable = "YearBuilt", min_lead = min(YearBuilt), max_lead = max(YearBuilt))
sum3 <- house |>
summarise(Variable = "MoSold", min_lead = min(MoSold), max_lead = max(MoSold))
sum4 <- house |>
summarise(Variable = "YrSold", min_lead = min(YrSold), max_lead = max(YrSold))
rbind(sum1, sum2, sum3, sum4)
Variable min_lead max_lead
1 YearRemodAdd 1950 2010
2 YearBuilt 1872 2010
3 MoSold 1 12
4 YrSold 2006 2010
Las cuatro variables miden instantes temporales (años o meses) sobre distintas cuestiones relacionadas con la vida de la vivienda. La variable MoSold toma valores del 1 al 12 porque informa únicamente del mes del año en el que se vendió la vivienda. Graficaremos a continuación las cuatro variables para comprobar como se distribuyen en función de nuestra objetivo SalePrice.
t1 <- house[!is.na(house$SalePrice), ] |>
ggplot(aes(x = YearRemodAdd, y = SalePrice)) +
geom_point(col = "#EB9891") +
geom_smooth(method = "lm", se = FALSE, color = "black", aes(group = 1)) +
scale_y_continuous(breaks= seq(0, 800000, by = 100000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
scale_x_continuous(labels = comma_format(big.mark = " ", decimal.mark = ",")) +
theme_minimal()
t2 <- house[!is.na(house$SalePrice), ] |>
ggplot(aes(x = YearBuilt, y = SalePrice)) +
geom_point(col = "#EB9891") +
geom_smooth(method = "lm", se = FALSE, color = "black", aes(group = 1)) +
scale_y_continuous(breaks= seq(0, 800000, by = 100000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
scale_x_continuous(labels = comma_format(big.mark = " ", decimal.mark = ",")) +
theme_minimal()
t3 <- aggregate(SalePrice ~ YrSold, house[!is.na(house$SalePrice), ], median) |>
mutate(n = pull(count(house[!is.na(house$SalePrice), ]$YrSold))) |>
ggplot(aes(x = YrSold, y = SalePrice)) +
geom_bar(stat = "identity", fill= "#EB9891") +
geom_label(aes(label = n, y = 0.1)) +
scale_y_continuous(breaks= seq(0, 800000, by = 25000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
coord_cartesian(ylim = c(0, 200000)) +
geom_hline(aes(yintercept = mean(SalePrice, na.rm = T), linetype = "Media"), colour = "black", size = .8) +
geom_hline(aes(yintercept = median(SalePrice, na.rm = T), linetype = "Mediana"), colour = "black", size = .8) +
scale_linetype_manual(name = "Medidas", values = c(Media = "solid", Mediana = "dotted")) +
labs(x = "Year Sold",y = "Sale Price") +
theme_minimal()
t4 <- aggregate(SalePrice ~ MoSold, house[!is.na(house$SalePrice), ], median) |>
mutate(n = pull(count(house[!is.na(house$SalePrice), ]$MoSold))) |>
ggplot(aes(x = MoSold, y = SalePrice)) +
geom_bar(stat = "identity", fill= "#EB9891") +
geom_label(aes(label = n, y = 0.1)) +
scale_y_continuous(breaks= seq(0, 800000, by = 25000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
coord_cartesian(ylim = c(0, 200000)) +
geom_hline(aes(yintercept = mean(SalePrice, na.rm = T), linetype = "Media"), colour = "black", size = .8) +
geom_hline(aes(yintercept = median(SalePrice, na.rm = T), linetype = "Mediana"), colour = "black", size = .8) +
scale_linetype_manual(name = "Medidas", values = c(Media = "solid", Mediana = "dotted")) +
labs(x = "Month Sold",y = "Sale Price") +
theme_minimal()
multiplot(t1, t2, t3, t4, cols = 2)
Algunos comentarios a la vista de los gráficos:
Para las variables YearRemodAdd y YearBuilt la correlación es claramente positiva: obviamente, a menos años pasen desde el momento de la reforma o construcción, mayor será el precio de venta de la vivienda. Además, a partir de estas dos variables podríamos construir otras que puedan aportar información también relevante a efectos de la modelización. Por un lado, la variable YearBuilt podría restarse con YrSold (\(YrSold-YearBuilt\)) para dar lugar a una nueva variable que ofrecería información acerca de la edad de la vivienda. Esta variable actuaría como penalizador a través de su correlación negativa con nuestra variable objetivo, pues a más edad, menor sería el precio de la vivienda. Además, la variable YearRemodAdd, al igual que como haremos con las variables relacionadas con las posibles piscinas de las viviendas, podría transformarse en binaria con las modalidades Yes y No y respondiendo a la pregunta de si se ha reformado o no en el tiempo.
Para las variables YrSold y MoSold no parece darse demasiada variación entre categorías. Para YrSold, los efectos de la crisis financiera de 2007 sí que parecen hacerse notar, sobre todo en la ligera disminución del precio de las viviendas a partir de ese año.
En resumen, con estas cuatro variables crearemos otras dos nuevas para, posteriormente, eliminar la variable YearRemodAdd original. Además, las variables YrSold y MoSold podrían funcionar también como cualitativas y factores.
sum1 <- house |>
summarise(Variable = "FullBath", min_lead = min(FullBath), max_lead = max(FullBath))
sum2 <- house |>
summarise(Variable = "HalfBath", min_lead = min(HalfBath), max_lead = max(HalfBath))
rbind(sum1, sum2)
Variable min_lead max_lead
1 FullBath 0 4
2 HalfBath 0 2
Estas dos variables cuentan los baños totales y los medios baños de la vivienda. Estas toman el valor 0 cuando la vivienda no dispone de estos tipos de baño. Graficaremos a continuación las dos variables para comprobar como se distribuyen en función de nuestra objetivo SalePrice.
c1 <- house[!is.na(house$SalePrice), ] |>
ggplot(aes(x = factor(FullBath), y = SalePrice)) +
geom_boxplot(col = "#EB9891") + labs(x = "Full Bath") +
scale_y_continuous(breaks= seq(0, 800000, by = 100000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
theme_minimal()
c2 <- house[!is.na(house$SalePrice), ] |>
ggplot(aes(x = factor(HalfBath), y = SalePrice)) +
geom_boxplot(col = "#EB9891") + labs(x = "Half Bath") +
scale_y_continuous(breaks= seq(0, 800000, by = 100000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
theme_minimal()
multiplot(c1, c2)
Como ninguna de estas dos variables parece muy potente, ni tampoco tiene una influencia determinante sobre nuestra objetivo, lo que se ha decidido es unirlas para conformar una nueva variable. De esta manera, la suma de FullBath y \(HalfBath * 0.5\) nos ofrecería información acerca del total de baños en la vivienda. La variable HalfBath (medios baños) se ponderaría por su significado real dividiendo sus valores y, por tanto, su importancia entre la mitad.
sum1 <- house |>
summarise(Variable = "BedroomAbvGr", min_lead = min(BedroomAbvGr), max_lead = max(BedroomAbvGr))
sum2 <- house |>
drop_na(GarageCars) |>
summarise(Variable = "GarageCars", min_lead = min(GarageCars), max_lead = max(GarageCars))
rbind(sum1, sum2)
Variable min_lead max_lead
1 BedroomAbvGr 0 8
2 GarageCars 0 5
Estas dos variables cuentan las habitaciones y los estacionamientos del garaje de la vivienda. Estas toman el valor 0 cuando la vivienda no dispone de estas áreas. Graficaremos a continuación las dos variables para comprobar como se distribuyen en función de nuestra objetivo SalePrice.
c1 <- house[!is.na(house$SalePrice), ] |>
ggplot(aes(x = factor(BedroomAbvGr), y = SalePrice)) +
geom_boxplot(col = "#EB9891") + labs(x = "Bedrooms") +
scale_y_continuous(breaks= seq(0, 800000, by = 100000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
theme_minimal()
c2 <- house[!is.na(house$SalePrice), ] |>
ggplot(aes(x = factor(GarageCars), y = SalePrice)) +
geom_boxplot(col = "#EB9891") + labs(x = "Parking lots") +
scale_y_continuous(breaks= seq(0, 800000, by = 100000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
theme_minimal()
multiplot(c1, c2)
Visualicemos también el peso de cada modalidad sobre el total de registros.
b1 <- ggplot(house, aes(x = factor(BedroomAbvGr))) +
geom_histogram(stat = "Count", fill = "#EB9891") +
labs(x="Bedroom") +
theme_minimal()
b2 <- ggplot(house, aes(x = factor(GarageCars))) +
geom_histogram(stat = "Count", fill = "#EB9891") +
labs(x="Parking lots") +
theme_minimal()
multiplot(b1, b2)
A la luz de los gráficos, se puede observar como ambas variables siguen una correlación positiva para con nuestra variable objetivo (conclusión coherente, como con las anteriores variables). A partir de los valores 3-4, el número de registros decae para esas modalidades. Seguramente se opte en la fase de reagrupación por eliminar o incluir esas modalidades en la categoría contigua.
A continuación analizaremos y agruparemos las variables cualitativas de cara a la creación de la receta del algoritmo KNN en modo regresión.
Rows: 2,919
Columns: 23
$ MSSubClass <fct> "2 story 1946+", "1 story 1946+", "2 story 194…
$ MSZoning <fct> RL, RL, RL, RL, RL, RL, RL, RL, RM, RL, RL, RL…
$ Street <fct> Pave, Pave, Pave, Pave, Pave, Pave, Pave, Pave…
$ Alley <fct> None, None, None, None, None, None, None, None…
$ LandContour <fct> Lvl, Lvl, Lvl, Lvl, Lvl, Lvl, Lvl, Lvl, Lvl, L…
$ Utilities <fct> AllPub, AllPub, AllPub, AllPub, AllPub, AllPub…
$ LandSlope <fct> Gtl, Gtl, Gtl, Gtl, Gtl, Gtl, Gtl, Gtl, Gtl, G…
$ Neighborhood <fct> CollgCr, Veenker, CollgCr, Crawfor, NoRidge, M…
$ BldgType <fct> 1Fam, 1Fam, 1Fam, 1Fam, 1Fam, 1Fam, 1Fam, 1Fam…
$ HouseStyle <fct> 2Story, 1Story, 2Story, 2Story, 2Story, 1.5Fin…
$ RoofStyle <fct> Gable, Gable, Gable, Gable, Gable, Gable, Gabl…
$ RoofMatl <fct> CompShg, CompShg, CompShg, CompShg, CompShg, C…
$ Exterior1st <fct> VinylSd, MetalSd, VinylSd, Wd Sdng, VinylSd, V…
$ Heating <fct> GasA, GasA, GasA, GasA, GasA, GasA, GasA, GasA…
$ CentralAir <fct> Y, Y, Y, Y, Y, Y, Y, Y, Y, Y, Y, Y, Y, Y, Y, Y…
$ Electrical <fct> SBrkr, SBrkr, SBrkr, SBrkr, SBrkr, SBrkr, SBrk…
$ Functional <fct> Typ, Typ, Typ, Typ, Typ, Typ, Typ, Typ, Min1, …
$ PavedDrive <fct> Y, Y, Y, Y, Y, Y, Y, Y, Y, Y, Y, Y, Y, Y, Y, Y…
$ Fence <fct> None, None, None, None, None, MnPrv, None, Non…
$ MoSold <fct> 2, 5, 9, 2, 12, 10, 8, 11, 4, 1, 2, 7, 9, 8, 5…
$ YrSold <fct> 2008, 2007, 2008, 2006, 2008, 2009, 2007, 2009…
$ SaleType <fct> WD, WD, WD, WD, WD, WD, WD, WD, WD, WD, WD, Ne…
$ SaleCondition <fct> Normal, Normal, Normal, Abnorml, Normal, Norma…
b1 <- aggregate(SalePrice ~ Street, house[!is.na(house$SalePrice), ], median) |>
mutate(n = pull(count(house[!is.na(house$SalePrice), ]$Street))) |>
ggplot(aes(x = Street, y = SalePrice)) +
geom_bar(stat = "identity", fill= "#56BCC2") +
geom_label(aes(label = n, y = 0.1)) +
scale_y_continuous(breaks= seq(0, 700000, by = 50000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_hline(aes(yintercept = mean(SalePrice, na.rm = T), linetype = "Media"), colour = "black", size = .8) +
geom_hline(aes(yintercept = median(SalePrice, na.rm = T), linetype = "Mediana"), colour = "black", size = .8) +
scale_linetype_manual(name = "Medidas", values = c(Media = "solid", Mediana = "dotted")) +
labs(x = "Street",y = "Sale Price") +
theme_minimal()
b2 <- aggregate(SalePrice ~ Alley, house[!is.na(house$SalePrice), ], median) |>
mutate(n = pull(count(house[!is.na(house$SalePrice), ]$Alley))) |>
ggplot(aes(x = Alley, y = SalePrice)) +
geom_bar(stat = "identity", fill= "#56BCC2") +
geom_label(aes(label = n, y = 0.1)) +
scale_y_continuous(breaks= seq(0, 700000, by = 50000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_hline(aes(yintercept = mean(SalePrice, na.rm = T), linetype = "Media"), colour = "black", size = .8) +
geom_hline(aes(yintercept = median(SalePrice, na.rm = T), linetype = "Mediana"), colour = "black", size = .8) +
scale_linetype_manual(name = "Medidas", values = c(Media = "solid", Mediana = "dotted")) +
labs(x = "Alley",y = "Sale Price") +
theme_minimal()
b3 <- aggregate(SalePrice ~ LandContour, house[!is.na(house$SalePrice), ], median) |>
mutate(n = pull(count(house[!is.na(house$SalePrice), ]$LandContour))) |>
ggplot(aes(x = LandContour, y = SalePrice)) +
geom_bar(stat = "identity", fill= "#56BCC2") +
geom_label(aes(label = n, y = 0.1)) +
scale_y_continuous(breaks= seq(0, 700000, by = 50000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_hline(aes(yintercept = mean(SalePrice, na.rm = T), linetype = "Media"), colour = "black", size = .8) +
geom_hline(aes(yintercept = median(SalePrice, na.rm = T), linetype = "Mediana"), colour = "black", size = .8) +
scale_linetype_manual(name = "Medidas", values = c(Media = "solid", Mediana = "dotted")) +
labs(x = "LandContour",y = "Sale Price") +
theme_minimal()
b4 <- aggregate(SalePrice ~ LandSlope, house[!is.na(house$SalePrice), ], median) |>
mutate(n = pull(count(house[!is.na(house$SalePrice), ]$LandSlope))) |>
ggplot(aes(x = LandSlope, y = SalePrice)) +
geom_bar(stat = "identity", fill= "#56BCC2") +
geom_label(aes(label = n, y = 0.1)) +
scale_y_continuous(breaks= seq(0, 700000, by = 50000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_hline(aes(yintercept = mean(SalePrice, na.rm = T), linetype = "Media"), colour = "black", size = .8) +
geom_hline(aes(yintercept = median(SalePrice, na.rm = T), linetype = "Mediana"), colour = "black", size = .8) +
scale_linetype_manual(name = "Medidas", values = c(Media = "solid", Mediana = "dotted")) +
labs(x = "LandSlope",y = "Sale Price") +
theme_minimal()
b5 <- aggregate(SalePrice ~ RoofStyle, house[!is.na(house$SalePrice), ], median) |>
mutate(n = pull(count(house[!is.na(house$SalePrice), ]$RoofStyle))) |>
ggplot(aes(x = RoofStyle, y = SalePrice)) +
geom_bar(stat = "identity", fill= "#56BCC2") +
geom_label(aes(label = n, y = 0.1)) +
scale_y_continuous(breaks= seq(0, 700000, by = 50000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_hline(aes(yintercept = mean(SalePrice, na.rm = T), linetype = "Media"), colour = "black", size = .8) +
geom_hline(aes(yintercept = median(SalePrice, na.rm = T), linetype = "Mediana"), colour = "black", size = .8) +
scale_linetype_manual(name = "Medidas", values = c(Media = "solid", Mediana = "dotted")) +
labs(x = "RoofStyle",y = "Sale Price") +
scale_x_discrete(guide = guide_axis(n.dodge = 2)) +
theme_minimal()
b6 <- aggregate(SalePrice ~ RoofMatl, house[!is.na(house$SalePrice), ], median) |>
mutate(n = pull(count(house[!is.na(house$SalePrice), ]$RoofMatl))) |>
ggplot(aes(x = RoofMatl, y = SalePrice)) +
geom_bar(stat = "identity", fill= "#56BCC2") +
geom_label(aes(label = n, y = 0.7)) +
scale_y_continuous(breaks= seq(0, 700000, by = 50000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_hline(aes(yintercept = mean(SalePrice, na.rm = T), linetype = "Media"), colour = "black", size = .8) +
geom_hline(aes(yintercept = median(SalePrice, na.rm = T), linetype = "Mediana"), colour = "black", size = .8) +
scale_linetype_manual(name = "Medidas", values = c(Media = "solid", Mediana = "dotted")) +
labs(x = "RoofMatl",y = "Sale Price") +
coord_flip() +
theme_minimal()
b7 <- aggregate(SalePrice ~ Exterior1st, house[!is.na(house$SalePrice), ], median) |>
mutate(n = pull(count(house[!is.na(house$SalePrice), ]$Exterior1st))) |>
ggplot(aes(x = Exterior1st, y = SalePrice)) +
geom_bar(stat = "identity", fill= "#56BCC2") +
geom_label(aes(label = n, y = 0.5)) +
scale_y_continuous(breaks= seq(0, 700000, by = 50000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_hline(aes(yintercept = mean(SalePrice, na.rm = T), linetype = "Media"), colour = "black", size = .8) +
geom_hline(aes(yintercept = median(SalePrice, na.rm = T), linetype = "Mediana"), colour = "black", size = .8) +
scale_linetype_manual(name = "Medidas", values = c(Media = "solid", Mediana = "dotted")) +
labs(x = "Exterior1st",y = "Sale Price") +
coord_flip() +
theme_minimal()
b8 <- aggregate(SalePrice ~ PavedDrive, house[!is.na(house$SalePrice), ], median) |>
mutate(n = pull(count(house[!is.na(house$SalePrice), ]$PavedDrive))) |>
ggplot(aes(x = PavedDrive, y = SalePrice)) +
geom_bar(stat = "identity", fill= "#56BCC2") +
geom_label(aes(label = n, y = 0.1)) +
scale_y_continuous(breaks= seq(0, 700000, by = 50000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_hline(aes(yintercept = mean(SalePrice, na.rm = T), linetype = "Media"), colour = "black", size = .8) +
geom_hline(aes(yintercept = median(SalePrice, na.rm = T), linetype = "Mediana"), colour = "black", size = .8) +
scale_linetype_manual(name = "Medidas", values = c(Media = "solid", Mediana = "dotted")) +
labs(x = "PavedDrive",y = "Sale Price") +
theme_minimal()
b9 <- aggregate(SalePrice ~ Fence, house[!is.na(house$SalePrice), ], median) |>
mutate(n = pull(count(house[!is.na(house$SalePrice), ]$Fence))) |>
ggplot(aes(x = Fence, y = SalePrice)) +
geom_bar(stat = "identity", fill= "#56BCC2") +
geom_label(aes(label = n, y = 0.1)) +
scale_y_continuous(breaks= seq(0, 700000, by = 50000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_hline(aes(yintercept = mean(SalePrice, na.rm = T), linetype = "Media"), colour = "black", size = .8) +
geom_hline(aes(yintercept = median(SalePrice, na.rm = T), linetype = "Mediana"), colour = "black", size = .8) +
scale_linetype_manual(name = "Medidas", values = c(Media = "solid", Mediana = "dotted")) +
labs(x = "Fence",y = "Sale Price") +
theme_minimal()
multiplot(b1, b2, b3, b4, b5, b6, b7, b8, b9, cols = 3)
A la luz de los gráficos, todas las variables presentan cierta influencia sobre nuestra variable objetivo, aunque muchas de sus categorías son muy minoritarias. Resulta necesario reagrupar la mayoría de variables. Se agruparan en función de su influencia sobre nuestra objetivo.
Street presenta dos categorías, pero tan solo una de ellas acapara la gran mayoría de registros (Pave). Se optará por eliminar esta variable.Alley presenta tres categorías, pero tan solo una de ellas acapara la gran mayoría de registros (None). Se optará por agrupar las otras dos Grvl y Pave en una sola.LandSlope presenta tres categorías, aunque dos tienen muy pocos registros. Se optará por agrupar las categorías Mod y Sev en una nueva.RoofStyle presenta seis categorías, aunque cuatro de ellas tienen registros inferiores a 10. Se optará por agrupar las categorías Flat, Shed y Mansard en la mayoritaria Hip. También se optará por agrupar la categoría Gambrel a Gable.RoofMatl presenta ocho categorías, aunque tan solo una tiene una proporción elevada de registros. Se optará por eliminar esta variable.LandContour presenta cuatro categorías, siendo una de ellas la mayoritaria. Se optará por agrupar las minoritarias Bnk, HLS y Low en una sola.Exterior1st presenta quince categorías. Agruparemos las categorías en función de la barrera de los 150 000 dólares por vivienda.PavedDrive presenta tres categorías, siendo una de ellas la mayoritaria. Se optará por agrupar las minoritarias N y P en una sola.Fence presenta cinco categorías. Se optará por agrupar las minoritarias GdWo, MnWw con la mayoritaria MnPrv en una sola. Por otro lado, también se optará por agrupar GdPrv con None.b1 <- aggregate(SalePrice ~ MSSubClass, house[!is.na(house$SalePrice), ], median) |>
mutate(n = pull(count(house[!is.na(house$SalePrice), ]$MSSubClass))) |>
ggplot(aes(x = MSSubClass, y = SalePrice)) +
geom_bar(stat = "identity", fill= "#56BCC2") +
geom_label(aes(label = n, y = 0.7)) +
scale_y_continuous(breaks= seq(0, 700000, by = 50000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_hline(aes(yintercept = mean(SalePrice, na.rm = T), linetype = "Media"), colour = "black", size = .8) +
geom_hline(aes(yintercept = median(SalePrice, na.rm = T), linetype = "Mediana"), colour = "black", size = .8) +
scale_linetype_manual(name = "Medidas", values = c(Media = "solid", Mediana = "dotted")) +
labs(x = "MSSubClass",y = "Sale Price") +
coord_flip() +
theme_minimal()
b2 <- aggregate(SalePrice ~ MSZoning, house[!is.na(house$SalePrice), ], median) |>
mutate(n = pull(count(house[!is.na(house$SalePrice), ]$MSZoning))) |>
ggplot(aes(x = MSZoning, y = SalePrice)) +
geom_bar(stat = "identity", fill= "#56BCC2") +
geom_label(aes(label = n, y = 0.1)) +
scale_y_continuous(breaks= seq(0, 700000, by = 50000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_hline(aes(yintercept = mean(SalePrice, na.rm = T), linetype = "Media"), colour = "black", size = .8) +
geom_hline(aes(yintercept = median(SalePrice, na.rm = T), linetype = "Mediana"), colour = "black", size = .8) +
scale_linetype_manual(name = "Medidas", values = c(Media = "solid", Mediana = "dotted")) +
labs(x = "MSZoning",y = "Sale Price") +
theme_minimal()
b3 <- aggregate(SalePrice ~ Neighborhood, house[!is.na(house$SalePrice), ], median) |>
mutate(n = pull(count(house[!is.na(house$SalePrice), ]$Neighborhood))) |>
ggplot(aes(x = Neighborhood, y = SalePrice)) +
geom_bar(stat = "identity", fill= "#56BCC2") +
geom_label(aes(label = n, y = 0.7)) +
scale_y_continuous(breaks= seq(0, 700000, by = 50000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_hline(aes(yintercept = mean(SalePrice, na.rm = T), linetype = "Media"), colour = "black", size = .8) +
geom_hline(aes(yintercept = median(SalePrice, na.rm = T), linetype = "Mediana"), colour = "black", size = .8) +
scale_linetype_manual(name = "Medidas", values = c(Media = "solid", Mediana = "dotted")) +
labs(x = "Neighborhood",y = "Sale Price") +
coord_flip() +
theme_minimal()
b4 <- aggregate(SalePrice ~ BldgType, house[!is.na(house$SalePrice), ], median) |>
mutate(n = pull(count(house[!is.na(house$SalePrice), ]$BldgType))) |>
ggplot(aes(x = BldgType, y = SalePrice)) +
geom_bar(stat = "identity", fill= "#56BCC2") +
geom_label(aes(label = n, y = 0.1)) +
scale_y_continuous(breaks= seq(0, 700000, by = 50000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_hline(aes(yintercept = mean(SalePrice, na.rm = T), linetype = "Media"), colour = "black", size = .8) +
geom_hline(aes(yintercept = median(SalePrice, na.rm = T), linetype = "Mediana"), colour = "black", size = .8) +
scale_linetype_manual(name = "Medidas", values = c(Media = "solid", Mediana = "dotted")) +
labs(x = "BldgType",y = "Sale Price") +
theme_minimal()
b5 <- aggregate(SalePrice ~ HouseStyle, house[!is.na(house$SalePrice), ], median) |>
mutate(n = pull(count(house[!is.na(house$SalePrice), ]$HouseStyle))) |>
ggplot(aes(x = HouseStyle, y = SalePrice)) +
geom_bar(stat = "identity", fill= "#56BCC2") +
geom_label(aes(label = n, y = 0.1)) +
scale_y_continuous(breaks= seq(0, 700000, by = 50000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_hline(aes(yintercept = mean(SalePrice, na.rm = T), linetype = "Media"), colour = "black", size = .8) +
geom_hline(aes(yintercept = median(SalePrice, na.rm = T), linetype = "Mediana"), colour = "black", size = .8) +
scale_linetype_manual(name = "Medidas", values = c(Media = "solid", Mediana = "dotted")) +
labs(x = "HouseStyle",y = "Sale Price") +
theme_minimal()
multiplot(b1, b2, b3, b4, b5, cols = 2)
MSSubClass presenta quince categorías muy distribuidas. Las agruparemos en tres categorías en función de su relación con nuestra variable objetivo. Distinguimos tres posibles subgrupos: \(< 150 000\), y \(> 150 000\).MSZoning presenta cinco categorías, pero tan solo dos de ellas acaparan la gran mayoría de registros (RL y RM). Se optará por agrupar las categorías minoritarias C (all) y RH en RM, y la categoría FV en RL.Neighborhood presenta veinticinco categorías muy distribuidas. Se optará por agrupar los barrios con viviendas más baratas y más caras juntos. Distinguimos cuatro posibles subgrupos: \(< 125 000\), \(125 000 > x > 150 000\), \(150 000 > x > 200 000\) y \(> 200 000\).BldgType presenta cinco categorías, pero tan solo dos de ellas acaparan la gran mayoría de registros (1Fam y TwnhsE). Se optará por agrupar las categorías mayoritarias 1Fam y TwnhsE juntas. Por su parte, las tres restantes también se agruparán en una nueva (2fmCon, Duplex y Twnhs).HouseStyle presenta ocho categorías, con tres siendo las mayoritarias en número de registros. Se optará por agrupar las categorías 2Story, 2.5Fin y SLvl en una nueva. Por su parte, también se creará una nueva modalidad para las categorías 1.5Fin, 1.5Unf, 1Story, 2.5Unf y SFoyer.b1 <- aggregate(SalePrice ~ Utilities, house[!is.na(house$SalePrice), ], median) |>
mutate(n = pull(count(house[!is.na(house$SalePrice), ]$Utilities))) |>
ggplot(aes(x = Utilities, y = SalePrice)) +
geom_bar(stat = "identity", fill= "#56BCC2") +
geom_label(aes(label = n, y = 0.1)) +
scale_y_continuous(breaks= seq(0, 700000, by = 50000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_hline(aes(yintercept = mean(SalePrice, na.rm = T), linetype = "Media"), colour = "black", size = .8) +
geom_hline(aes(yintercept = median(SalePrice, na.rm = T), linetype = "Mediana"), colour = "black", size = .8) +
scale_linetype_manual(name = "Medidas", values = c(Media = "solid", Mediana = "dotted")) +
labs(x = "Utilities",y = "Sale Price") +
theme_minimal()
b2 <- aggregate(SalePrice ~ Heating, house[!is.na(house$SalePrice), ], median) |>
mutate(n = pull(count(house[!is.na(house$SalePrice), ]$Heating))) |>
ggplot(aes(x = Heating, y = SalePrice)) +
geom_bar(stat = "identity", fill= "#56BCC2") +
geom_label(aes(label = n, y = 0.1)) +
scale_y_continuous(breaks= seq(0, 700000, by = 50000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_hline(aes(yintercept = mean(SalePrice, na.rm = T), linetype = "Media"), colour = "black", size = .8) +
geom_hline(aes(yintercept = median(SalePrice, na.rm = T), linetype = "Mediana"), colour = "black", size = .8) +
scale_linetype_manual(name = "Medidas", values = c(Media = "solid", Mediana = "dotted")) +
labs(x = "Heating",y = "Sale Price") +
theme_minimal()
b3 <- aggregate(SalePrice ~ CentralAir, house[!is.na(house$SalePrice), ], median) |>
mutate(n = pull(count(house[!is.na(house$SalePrice), ]$CentralAir))) |>
ggplot(aes(x = CentralAir, y = SalePrice)) +
geom_bar(stat = "identity", fill= "#56BCC2") +
geom_label(aes(label = n, y = 0.1)) +
scale_y_continuous(breaks= seq(0, 700000, by = 50000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_hline(aes(yintercept = mean(SalePrice, na.rm = T), linetype = "Media"), colour = "black", size = .8) +
geom_hline(aes(yintercept = median(SalePrice, na.rm = T), linetype = "Mediana"), colour = "black", size = .8) +
scale_linetype_manual(name = "Medidas", values = c(Media = "solid", Mediana = "dotted")) +
labs(x = "CentralAir",y = "Sale Price") +
theme_minimal()
b4 <- aggregate(SalePrice ~ Electrical, house[!is.na(house$SalePrice), ], median) |>
ggplot(aes(x = Electrical, y = SalePrice)) +
geom_bar(stat = "identity", fill= "#56BCC2") +
scale_y_continuous(breaks= seq(0, 700000, by = 50000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_hline(aes(yintercept = mean(SalePrice, na.rm = T), linetype = "Media"), colour = "black", size = .8) +
geom_hline(aes(yintercept = median(SalePrice, na.rm = T), linetype = "Mediana"), colour = "black", size = .8) +
scale_linetype_manual(name = "Medidas", values = c(Media = "solid", Mediana = "dotted")) +
labs(x = "Electrical",y = "Sale Price") +
theme_minimal()
b5 <- aggregate(SalePrice ~ Functional, house[!is.na(house$SalePrice), ], median) |>
mutate(n = pull(count(house[!is.na(house$SalePrice), ]$Functional))) |>
ggplot(aes(x = Functional, y = SalePrice)) +
geom_bar(stat = "identity", fill= "#56BCC2") +
geom_label(aes(label = n, y = 0.1)) +
scale_y_continuous(breaks= seq(0, 700000, by = 50000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_hline(aes(yintercept = mean(SalePrice, na.rm = T), linetype = "Media"), colour = "black", size = .8) +
geom_hline(aes(yintercept = median(SalePrice, na.rm = T), linetype = "Mediana"), colour = "black", size = .8) +
scale_linetype_manual(name = "Medidas", values = c(Media = "solid", Mediana = "dotted")) +
labs(x = "Functional",y = "Sale Price") +
theme_minimal()
multiplot(b1, b2, b3, b4, b5, cols = 2)
Utilities presenta dos categorías, pero tan solo una de ellas acapara todos los registros salvo uno. Se optará por eliminar esta variable.Heating presenta seis categorías, pero tan solo dos de ellas acaparan la gran mayoría de registros. Se optará por eliminar esta variable.CentralAir presenta tan solo dos categorías. La mantendremos tal y como está.Electrical presenta cinco categorías. En este caso, se optará por agrupar todas las minoritarias en una nueva y mantener la mayoritaria SBrkr.Functional presenta siete categorías, con una de ellas siendo la mayoritaria en número de registros. Se optará por agrupar todas las minoritarias en una nueva y mantener la mayoritaria Typ.b1 <- aggregate(SalePrice ~ SaleType, house[!is.na(house$SalePrice), ], median) |>
mutate(n = pull(count(house[!is.na(house$SalePrice), ]$SaleType))) |>
ggplot(aes(x = SaleType, y = SalePrice)) +
geom_bar(stat = "identity", fill= "#56BCC2") +
geom_label(aes(label = n, y = 0.1)) +
scale_y_continuous(breaks= seq(0, 700000, by = 50000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_hline(aes(yintercept = mean(SalePrice, na.rm = T), linetype = "Media"), colour = "black", size = .8) +
geom_hline(aes(yintercept = median(SalePrice, na.rm = T), linetype = "Mediana"), colour = "black", size = .8) +
scale_linetype_manual(name = "Medidas", values = c(Media = "solid", Mediana = "dotted")) +
labs(x = "SaleType",y = "Sale Price") +
theme_minimal()
b2 <- aggregate(SalePrice ~ SaleCondition, house[!is.na(house$SalePrice), ], median) |>
mutate(n = pull(count(house[!is.na(house$SalePrice), ]$SaleCondition))) |>
ggplot(aes(x = SaleCondition, y = SalePrice)) +
geom_bar(stat = "identity", fill= "#56BCC2") +
geom_label(aes(label = n, y = 0.1)) +
scale_y_continuous(breaks= seq(0, 700000, by = 50000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_hline(aes(yintercept = mean(SalePrice, na.rm = T), linetype = "Media"), colour = "black", size = .8) +
geom_hline(aes(yintercept = median(SalePrice, na.rm = T), linetype = "Mediana"), colour = "black", size = .8) +
scale_linetype_manual(name = "Medidas", values = c(Media = "solid", Mediana = "dotted")) +
labs(x = "SaleCondition",y = "Sale Price") +
theme_minimal()
multiplot(b1, b2)
SaleType presenta nueve categorías, pero tan solo dos de ellas acaparan la gran mayoría de registros (WD y New). Se optará por agrupar las categorías con valores más altos respecto de nuestra objetivo, y viceversa. De esta manera, las categorías Con, CWD y New irán juntas. Las cinco minoritarias restantes se añadirán a la mayoritaria WD.SaleCondition presenta seis categorías, con tres siendo las mayoritarias en número de registros. Se optará por agrupar las categorías Normal, Abnorml, AdjLand, Alloca y Family en una nueva. Por otro lado, se mantendrá a la categoría Partial tal y como está.Tras la fase de exploración de los datos, continuaremos con las fases de muestreo y modificación de los datos. Dado que nuestro dataset contiene tan solo 2919 observaciones (1460 para train y 1459 para test), no será necesario realizar muestreo (nos quedaríamos prácticamente sin filas si lo hacemos). Por otro lado, iremos fijando semilla a fin de poder sacar conclusiones de los resultados de nuestros distintos modelos.
En segundo lugar, para la fase de modificación de los datos, consideraremos dos apartados principales. Uno primero en donde se ejecutarán las modificaciones estructurales que afecten a toda las base de datos (transformar variables a factores, problemas de codificación o de rango, variables que no aportan, creación de variables en general, etc.), y uno segundo en donde se llevarán a cabo aquellas modificaciones que afecten a cada algoritmo en concreto a modo de receta (normalización para la métrica, recategorización, tratamiento de outliers/ausentes, dummyficación, etc.).
En este caso, como las particiones ya están realizadas, vamos a aplicar al conjunto completo (train+ test) las modificaciones estructurales para que tengan la misma forma.
# Aplicamos las modificaciones estructurales al dataset completo (train + test)
## Tratamiento de ausentes y transformación de variables texto a factor
house_complete <-
house_complete |>
mutate(Alley =
ifelse(is.na(Alley), "None", Alley)) |>
mutate(BsmtQual =
ifelse(is.na(BsmtQual), "None", BsmtQual)) |>
mutate(BsmtCond =
ifelse(is.na(BsmtCond), "None", BsmtCond)) |>
mutate(FireplaceQu =
ifelse(is.na(FireplaceQu), "None", FireplaceQu)) |>
mutate(GarageQual =
ifelse(is.na(GarageQual), "None", GarageQual)) |>
mutate(PoolQC =
ifelse(is.na(PoolQC), "None", PoolQC)) |>
mutate(Fence =
ifelse(is.na(Fence), "None", Fence))
## Modificaciones en la tipología de las variables
house_complete$ExterCond <-
as.integer(revalue(house_complete$ExterCond, quality))
house_complete$BsmtQual <-
as.integer(revalue(house_complete$BsmtQual, quality))
house_complete$BsmtCond <-
as.integer(revalue(house_complete$BsmtCond, quality))
house_complete$HeatingQC <-
as.integer(revalue(house_complete$HeatingQC, quality))
house_complete$KitchenQual <-
as.integer(revalue(house_complete$KitchenQual, quality))
house_complete$FireplaceQu <-
as.integer(revalue(house_complete$FireplaceQu, quality))
house_complete$GarageQual <-
as.integer(revalue(house_complete$GarageQual, quality))
house_complete$PoolQC <-
as.integer(revalue(house_complete$PoolQC, quality))
house_complete$KitchenQual[is.na(house_complete$KitchenQual)] <- 3
house_complete$GarageCars[is.na(house_complete$GarageCars)] <- 2
## Creación de nuevas variables
house_complete$RemodAdd <- # variable binaria ¿reformada?
ifelse(house_complete$YearBuilt == house_complete$YearRemodAdd, 0, 1)
house_complete$New <- # variable binaria ¿es nueva?
ifelse(house_complete$YrSold == house_complete$YearBuilt, 1, 0)
house_complete$TotBath <- # variable baño total
house_complete$FullBath + (house_complete$HalfBath * 0.5)
house_complete <-
house_complete |>
select(-c(FullBath, HalfBath))
house_complete$Age <- # variable edad de la vivienda
house_complete$YrSold-house_complete$YearRemodAdd
house_complete <-
house_complete |>
select(-c(YearRemodAdd))
## Transformación del resto a factor
house_complete$MSSubClass <-
as.factor(house_complete$MSSubClass)
house_complete$YrSold <-
as.factor(house_complete$YrSold)
house_complete$MoSold <-
as.factor(house_complete$MoSold)
house_complete <-
house_complete |>
mutate_if(~!is.numeric(.), as.factor)
## Eliminación de variables por problemas de colinealidad
house_complete <-
house_complete |>
select(-c(GarageYrBlt, TotalBsmtSF, TotRmsAbvGrd, Fireplaces, GarageArea, PoolArea))
## Eliminación de variables por insustanciales
house_complete <-
house_complete |>
select(-c(Id, OpenPorchSF, Street, RoofMatl, Utilities, Heating, BedroomAbvGr, PoolQC))
Tras homogeneizar los subconjuntos de train y test, los volvemos a subdividir tal y como se encontraban en origen.
house_train <- house_complete[1:1460,]
house_test <- house_complete[1461:2919,]
Tras las particiones, definimos la receta, indicándole el conjunto donde tenemos validación y train, y enfrentamos nuestra variable objetivo SalePrice con todas las demás.
Después, asignamos posibles roles, sujetos a modificación, que nos permitan diferenciar acciones entre las variables.
# Receta
house_rec_knn <-
# Fórmula y datos
recipe(data = house_train, SalePrice ~ .)|>
# Roles
add_role(where(is.factor),
new_role = "cualitativa") |>
add_role(where(is.numeric),
new_role = "cuantitativa") |>
add_role(RemodAdd, New,
new_role = "binaria") |>
add_role(LotFrontage, LotArea, `1stFlrSF`, `2ndFlrSF`, GrLivArea, Age,
new_role = "mediana") |>
add_role(all_nominal_predictors(),
new_role = "moda") |>
add_role(ExterCond, BsmtCond, HeatingQC, GarageQual, KitchenQual,
new_role = "imputar mediana") |>
add_role(OverallCond, BsmtQual, FireplaceQu, TotBath, YearBuilt, GarageCars,
new_role = "imputar media")
En este apartado, reagrupamos las variables cualitativas con step_mutate con lo definido en la fase de exploración. Para el algoritmo de clasificación KNN, cuanto más agrupadas estén este tipo de categorías, mejor.
house_rec_knn <-
house_rec_knn |>
step_mutate(Alley =
fct_collapse(Alley,
Grvl_Pave = c("Grvl", "Pave"))) |>
step_mutate(LandSlope =
fct_collapse(LandSlope,
Mod_Sev = c("Mod", "Sev"))) |>
step_mutate(LandContour =
fct_collapse(LandContour,
Bnk_HLS_Low = c("Bnk", "HLS", "Low"))) |>
step_mutate(RoofStyle =
fct_collapse(RoofStyle,
Hip = c("Flat", "Shed", "Mansard", "Hip"),
Gam_Gable = c("Gambrel", "Gable"))) |>
step_mutate(Exterior1st =
fct_collapse(Exterior1st,
`<150k` = c("AsphShn", "BrkComm", "BrkFace",
"CemntBd", "HdBoard", "MetalSd",
"Stucco", "Wd Sdng", "WdShing",
"AsbShng"),
`>150k` = c("CBlock", "ImStucc", "Plywood",
"Stone", "VinylSd"))) |>
step_mutate(PavedDrive =
fct_collapse(PavedDrive,
N_P = c("N", "P"))) |>
step_mutate(Fence =
fct_collapse(Fence,
MnPrv = c("GdWo", "MnWw"),
None = c("GdPrv"))) |>
step_mutate(MSSubClass =
fct_collapse(MSSubClass,
`<150k` = c("30", "40", "45", "50",
"80", "85", "120", "160",
"180", "190"),
`>150k` = c("20", "60", "70", "75",
"90", "150"))) |>
step_mutate(MSZoning =
fct_collapse(MSZoning,
RM = c("C (all)", "RH"),
RL = c("FV"))) |>
step_mutate(Neighborhood =
fct_collapse(Neighborhood,
`<125k` = c("BrDale", "BrkSide", "Edwards",
"IDOTRR", "MeadowV", "OldTown"),
`125k-150k` = c("Blmngtn", "NAmes", "NPkVill",
"SWISU", "Sawyer"),
`150k-200k` = c("Blueste", "CollgCr", "Gilbert",
"Mitchel", "NWAmes", "SawyerW"),
`>200k` = c("ClearCr", "Crawfor", "NoRidge",
"NridgHt", "Somerst", "StoneBr",
"Timber", "Veenker"))) |>
step_mutate(BldgType =
fct_collapse(BldgType,
Fam_TwnhsE = c("1Fam", "TwnhsE"),
fmCon_Dup_Twnhs = c("2fmCon", "Duplex", "Twnhs"))) |>
step_mutate(HouseStyle =
fct_collapse(HouseStyle,
Story_Fin_SLvl = c("2Story", "2.5Fin", "SLvl"),
fmCon_Dup_Twnhs = c("1.5Fin", "1.5Unf", "1Story",
"2.5Unf", "SFoyer"))) |>
step_mutate(Electrical =
fct_collapse(Electrical,
Other = c("FuseA", "FuseF", "FuseP", "Mix"))) |>
step_mutate(Functional =
fct_collapse(Functional,
Other = c("Maj1", "Maj2", "Min1", "Min2",
"Mod", "Sev"))) |>
step_mutate(SaleType =
fct_collapse(SaleType,
Con_CWD_New = c("Con", "CWD", "New"),
WD = c("COD", "ConLD", "ConLI", "ConLw",
"Oth"))) |>
step_mutate(SaleCondition =
fct_collapse(SaleCondition,
Normal_Others = c("Abnorml", "AdjLand", "Alloca",
"Family"))) |>
step_mutate(MoSold =
fct_collapse(MoSold,
Jan_Ap_May_Oct = c("1", "4", "5", "10"),
Mar_Jun_Jul = c("3", "6", "7"),
Feb_Aug_Sep_Nov_Dec = c("2", "8", "9",
"11", "12"))) |>
step_mutate(YrSold =
fct_collapse(YrSold,
`2006-2009` = c("2006", "2007", "2008", "2009")))
En el algoritmo de clasificación KNN, la recategorización de las variables cuantitativas no es estrictamente necesaria. Por tanto, en este receta, del total de variables cuantitativas no vamos a modificar ninguna.
box1 <-
ggplot(house_complete, aes(SalePrice, Age)) +
geom_boxplot() +
theme_minimal() +
coord_flip()
box2 <-
ggplot(house_complete, aes(SalePrice, LotFrontage)) +
geom_boxplot() +
theme_minimal() +
coord_flip()
box3 <-
ggplot(house_complete, aes(SalePrice, LotArea)) +
geom_boxplot() +
theme_minimal() +
coord_flip()
box4 <-
ggplot(house_complete, aes(SalePrice, OverallCond)) +
geom_boxplot() +
theme_minimal() +
coord_flip()
box5 <-
ggplot(house_complete, aes(SalePrice, YearBuilt)) +
geom_boxplot() +
theme_minimal() +
coord_flip()
box6 <-
ggplot(house_complete, aes(SalePrice, ExterCond)) +
geom_boxplot() +
theme_minimal() +
coord_flip()
box7 <-
ggplot(house_complete, aes(SalePrice, BsmtQual)) +
geom_boxplot() +
theme_minimal() +
coord_flip()
box8 <-
ggplot(house_complete, aes(SalePrice, BsmtCond)) +
geom_boxplot() +
theme_minimal() +
coord_flip()
box9 <-
ggplot(house_complete, aes(SalePrice, HeatingQC)) +
geom_boxplot() +
theme_minimal() +
coord_flip()
box10 <-
ggplot(house_complete, aes(SalePrice, `1stFlrSF`)) +
geom_boxplot() +
theme_minimal() +
coord_flip()
box11 <-
ggplot(house_complete, aes(SalePrice, `2ndFlrSF`)) +
geom_boxplot() +
theme_minimal() +
coord_flip()
box12 <-
ggplot(house_complete, aes(SalePrice, GrLivArea)) +
geom_boxplot() +
theme_minimal() +
coord_flip()
box13 <-
ggplot(house_complete, aes(SalePrice, KitchenQual)) +
geom_boxplot() +
theme_minimal() +
coord_flip()
box14 <-
ggplot(house_complete, aes(SalePrice, FireplaceQu)) +
geom_boxplot() +
theme_minimal() +
coord_flip()
box15 <-
ggplot(house_complete, aes(SalePrice, GarageCars)) +
geom_boxplot() +
theme_minimal() +
coord_flip()
box16 <-
ggplot(house_complete, aes(SalePrice, GarageQual)) +
geom_boxplot() +
theme_minimal() +
coord_flip()
box17 <-
ggplot(house_complete, aes(SalePrice, TotBath)) +
geom_boxplot() +
theme_minimal() +
coord_flip()
multiplot(box1, box2, box3, box4, box5, box6, box7, box8, box9, box10, cols = 2)
multiplot(box11, box12, box13, box14, box15, box16, box17, cols = 2)
Si observamos estos gráficos de cajas y bigotes, todas nuestras variables cuantitativas continuas son asimétricas, por lo que se detectarán los outliers y se imputarán los ausentes por la mediana. Para el caso de las semi-cuantitativas (KitchenQual, por ejemplo) se imputarán directamente los ausentes por uno de los dos métodos en función de la simetría o asimetría de la variable. Para el resto de variables cualitativas, imputamos directamente por la moda.
house_rec_knn <-
house_rec_knn |>
# Detección de outliers por la mediana y por la media
step_mutate(across(has_role("mediana"), function(x) { ifelse(abs(scores(x, type = "mad")) > 3 & !is.na(x), NA, x) })) |>
# Imputación de ausentes por la mediana, la media y la moda
step_impute_median(has_role("mediana")) |>
step_impute_median(has_role("imputar_mediana")) |>
step_impute_mean(has_role("imputar_media")) |>
step_impute_mode(has_role("moda"))
Aplicamos el filtro de correlación a nuestras variables numéricas.
En el caso del algoritmo KNN, necesitaremos normalizar nuestras variables por rango para que todas tengan el mismo peso, entre 0 y 1.
house_rec_knn <-
house_rec_knn |>
step_range(all_numeric_predictors(), min = 0, max = 1) |>
step_other(all_nominal_predictors(), threshold = 0.02)
Como esta receta es para el modelo KNN, debemos dummyficar nuestras variables cualitativas. Para ello, tomamos todas las nominales, menos nuestra variable objetivo.
house_rec_knn <-
house_rec_knn |>
step_dummy(all_nominal(), -all_outcomes())
house_rec_knn <-
house_rec_knn |>
step_zv(all_predictors())
Por último, horneamos nuestra receta para comprobar que todas nuestras nuevas variables recategorizadas se hayan creado correctamente.
# A tibble: 1,460 × 43
LotFrontage LotArea OverallCond YearBuilt ExterCond BsmtQual
<dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
1 0.151 0.0334 0.5 0.949 0.5 0.8
2 0.202 0.0388 0.875 0.754 0.5 0.8
3 0.161 0.0465 0.5 0.935 0.5 0.8
4 0.134 0.0386 0.5 0.312 0.5 0.6
5 0.216 0.0606 0.5 0.928 0.5 0.8
6 0.219 0.0599 0.5 0.877 0.5 0.8
7 0.185 0.0411 0.5 0.957 0.5 1
8 0.164 0.0425 0.625 0.732 0.5 0.8
9 0.103 0.0225 0.5 0.428 0.5 0.6
10 0.0993 0.0286 0.625 0.486 0.5 0.6
# … with 1,450 more rows, and 37 more variables: BsmtCond <dbl>,
# HeatingQC <dbl>, `1stFlrSF` <dbl>, `2ndFlrSF` <dbl>,
# GrLivArea <dbl>, KitchenQual <dbl>, FireplaceQu <dbl>,
# GarageCars <dbl>, GarageQual <dbl>, RemodAdd <dbl>, New <dbl>,
# TotBath <dbl>, Age <dbl>, SalePrice <dbl>,
# MSSubClass_X.150k <dbl>, MSZoning_RL <dbl>, Alley_None <dbl>,
# LandContour_Lvl <dbl>, LandSlope_Mod_Sev <dbl>, …
En este segundo apartado repetiremos la receta, pero esta vez adaptada a los árboles de clasificación.
# Receta
house_rec_arbol <-
# Fórmula y datos
recipe(data = house_train, SalePrice ~ .)|>
# Roles
add_role(where(is.factor), new_role = "cualitativa") |>
add_role(where(is.numeric), new_role = "cuantitativa") |>
add_role(RemodAdd, New, new_role = "binaria") |>
add_role(all_nominal_predictors(), new_role = "moda")
En este apartado, reagrupamos las variables cualitativas con step_mutate. Para el método de árboles de clasificación, no es necesario agrupar en exceso este tipo de variables, por lo que únicamente agruparemos las variables Neighborhood y Exterior1st (las que más modalidades presentaban).
house_rec_arbol <-
house_rec_arbol |>
step_mutate(Neighborhood =
fct_collapse(Neighborhood,
`<125k` = c("BrDale", "BrkSide", "Edwards",
"IDOTRR", "MeadowV", "OldTown"),
`125k-150k` = c("Blmngtn", "NAmes", "NPkVill",
"SWISU", "Sawyer"),
`150k-200k` = c("Blueste", "CollgCr", "Gilbert",
"Mitchel", "NWAmes", "SawyerW"),
`>200k` = c("ClearCr", "Crawfor", "NoRidge",
"NridgHt", "Somerst", "StoneBr",
"Timber", "Veenker"))) |>
step_mutate(Exterior1st =
fct_collapse(Exterior1st,
`<150k` = c("AsphShn", "BrkComm", "BrkFace",
"CemntBd", "HdBoard", "MetalSd",
"Stucco", "Wd Sdng", "WdShing",
"AsbShng"),
`>150k` = c("CBlock", "ImStucc", "Plywood",
"Stone", "VinylSd")))
Como esta receta es para el método de árboles de clasificación, tenemos que recategorizar todas las variables cuantitativas.
house_rec_arbol <-
house_rec_arbol |>
step_mutate(GrLivArea =
cut(GrLivArea,
breaks = c(0, 1000, 2000, 3000, Inf),
labels = c("<1000", "1000-2000", "2000-3000", ">3000"))) |>
step_mutate(`1stFlrSF` =
cut(`1stFlrSF`,
breaks = c(0, 1000, 2000, Inf),
labels = c("<1000", "1000-2000", ">2000"))) |>
step_mutate(`2ndFlrSF` =
cut(`2ndFlrSF`,
breaks = c(0, 500, 1000, 1500, Inf),
labels = c("<500", "500-1000", "1000-1500", ">1500"))) |>
step_mutate(LotArea =
cut(LotArea,
breaks = c(0, 10000, 20000, 30000, Inf),
labels = c("<10000", "10000-20000", "20000-30000", ">30000"))) |>
step_mutate(LotFrontage =
cut(LotFrontage,
breaks = c(0, 50, 100, 150, Inf),
labels = c("<50", "50-100", "100-150", ">150"))) |>
step_mutate(OverallCond =
cut(OverallCond,
breaks = c(-Inf, 3, 4, 5, 6, 7, Inf),
labels = c("1-3", "4", "5", "6", "7", "8-9"))) |>
step_mutate(BsmtCond =
cut(BsmtCond,
breaks = c(-Inf, 1, 2, 3, Inf),
labels = c("0-1", "2", "3", "4"))) |>
step_mutate(FireplaceQu =
cut(FireplaceQu,
breaks = c(-Inf, 2, 4, Inf),
labels = c("0-2", "3-4", "5"))) |>
step_mutate(ExterCond =
cut(ExterCond,
breaks = c(-Inf, 2, 3, Inf),
labels = c("0-2", "3", "4-5"))) |>
step_mutate(KitchenQual =
cut(KitchenQual,
breaks = c(-Inf, 3, 4, Inf),
labels = c("2-3", "4", "5"))) |>
step_mutate(BsmtQual =
cut(BsmtQual,
breaks = c(-Inf, 2, 3, 4, Inf),
labels = c("0&2", "3", "4", "5"))) |>
step_mutate(GarageQual =
cut(GarageQual,
breaks = c(-Inf, 1, 2, Inf),
labels = c("0-1", "2", "3-5"))) |>
step_mutate(HeatingQC =
cut(HeatingQC,
breaks = c(-Inf, 3, 4, Inf),
labels = c("1-3", "4", "5"))) |>
step_mutate(YearBuilt =
cut(YearBuilt,
breaks = c(-Inf, 1920, 1960, Inf),
labels = c("<1920", "1920-1960", ">1960"))) |>
step_mutate(GarageCars =
cut(GarageCars,
breaks = c(-Inf, 1, 2, 3, Inf),
labels = c("0-1", "2", "3", "4"))) |>
step_mutate(Age =
cut(Age,
breaks = c(-Inf, 10, 20, 30, 40, 50, Inf),
labels = c("0-10", "10-20", "20-30", "30-40", "40-50", "+50"))) |>
step_mutate(TotBath =
cut(TotBath,
breaks = c(-Inf, 1, 1.5, 2, 2.5, Inf),
labels = c("0-1", "1.5", "2", "2.5", "+3"))) |>
step_mutate(RemodAdd =
cut(RemodAdd,
breaks = c(-Inf, 0, Inf),
labels = c("No", "Yes"))) |>
step_mutate(New =
cut(New,
breaks = c(-Inf, 0, Inf),
labels = c("No", "Yes"))) |>
step_mutate(across(where(is.character), as_factor))
multiplot(box1, box2, box3, box4, box5, box6, box7, box8, box9, box10, cols = 2)
multiplot(box11, box12, box13, box14, box15, box16, box17, cols = 2)
En esta receta se ha optado por recategorizar y reagrupar todas las variables, incluso las cuantitativas. De esta manera, no hay outliers que detectar, pues todas las nuevas categorías se han creado manualmente a partir de la agrupación de las antiguas modalidades. Aún así, imputamos a los posibles ausentes la moda al haber convertido todas las variables en cualitativas.
house_rec_arbol <-
house_rec_arbol |>
step_impute_mode(all_predictors())
Como todas nuestras variables están ahora recategorizadas y convertidas a factor, no existen variables estrictamente numéricas a las que aplicar el filtro de correlación.
En árboles, no será necesario normalizar por rango.
En árboles, tampoco será necesario dummyficar.
house_rec_arbol <-
house_rec_arbol |>
step_zv(all_predictors())
Por último, horneamos nuestra receta para comprobar que todas nuestras nuevas variables recategorizadas se hayan creado correctamente.
# A tibble: 1,460 × 39
MSSubClass MSZoning LotFrontage LotArea Alley LandContour LandSlope
<fct> <fct> <fct> <fct> <fct> <fct> <fct>
1 60 RL 50-100 <10000 None Lvl Gtl
2 20 RL 50-100 <10000 None Lvl Gtl
3 60 RL 50-100 10000-… None Lvl Gtl
4 70 RL 50-100 <10000 None Lvl Gtl
5 60 RL 50-100 10000-… None Lvl Gtl
6 50 RL 50-100 10000-… None Lvl Gtl
7 20 RL 50-100 10000-… None Lvl Gtl
8 60 RL 50-100 10000-… None Lvl Gtl
9 50 RM 50-100 <10000 None Lvl Gtl
10 190 RL <50 <10000 None Lvl Gtl
# … with 1,450 more rows, and 32 more variables: Neighborhood <fct>,
# BldgType <fct>, HouseStyle <fct>, OverallCond <fct>,
# YearBuilt <fct>, RoofStyle <fct>, Exterior1st <fct>,
# ExterCond <fct>, BsmtQual <fct>, BsmtCond <fct>, HeatingQC <fct>,
# CentralAir <fct>, Electrical <fct>, `1stFlrSF` <fct>,
# `2ndFlrSF` <fct>, GrLivArea <fct>, KitchenQual <fct>,
# Functional <fct>, FireplaceQu <fct>, GarageCars <fct>, …
Una vez definida la receta, definimos el modelo y lo unimos con nuestra receta creando un flujo de clasificación.
Se trata de un nuevo modelo con tune() para poder definir posteriormente los parámetros de forma manual. Activamos, además, el modo regresión con mode = "regression"
# Modelo tuneado
knn_model_tune <-
nearest_neighbor(mode = "regression", neighbors = tune("k"),
weight_func = tune("weight"), dist_power = tune("dist")) |>
set_engine("kknn")
# Flujo de trabajo
house_wflow_knn <-
workflow() |>
add_recipe(house_rec_knn) |>
add_model(knn_model_tune)
Para esta práctica evaluaremos directamente a través del método de Validación Cruzada V-Folds, que nos proporcionará además una métrica media con la desviación típica para cada modelo.
Para ello, en primer lugar crearemos un grid con la función expand_grid que nos permitirá definir manualmente los parámetros que queramos sin estar condicionados por una función automática.
De esta manera, podremos probar todas las combinaciones que queramos.
Tras comprobar varias combinaciones, finalmente nos decantamos por la siguiente.
grid_knn <-
expand_grid("k" = c(8, 20, 60, 80, 100),
"weight" = c("inv", "gaussian"),
"dist" = c(0.01, 0.9, 2, 5, 15, 20))
grid_knn
# A tibble: 60 × 3
k weight dist
<dbl> <chr> <dbl>
1 8 inv 0.01
2 8 inv 0.9
3 8 inv 2
4 8 inv 5
5 8 inv 15
6 8 inv 20
7 8 gaussian 0.01
8 8 gaussian 0.9
9 8 gaussian 2
10 8 gaussian 5
# … with 50 more rows
Aplicamos en esta fase el proceso de Validación Cruzada V-Folds.
house_cv_folds <- vfold_cv(data = house_train, v = 8, repeats = 4)
house_cv_folds
# 8-fold cross-validation repeated 4 times
# A tibble: 32 × 3
splits id id2
<list> <chr> <chr>
1 <split [1277/183]> Repeat1 Fold1
2 <split [1277/183]> Repeat1 Fold2
3 <split [1277/183]> Repeat1 Fold3
4 <split [1277/183]> Repeat1 Fold4
5 <split [1278/182]> Repeat1 Fold5
6 <split [1278/182]> Repeat1 Fold6
7 <split [1278/182]> Repeat1 Fold7
8 <split [1278/182]> Repeat1 Fold8
9 <split [1277/183]> Repeat2 Fold1
10 <split [1277/183]> Repeat2 Fold2
# … with 22 more rows
set.seed(05492)
clusters <- detectCores() - 1
make_cluster <- makeCluster(clusters)
registerDoParallel(make_cluster)
house_knn_vf <-
house_wflow_knn |>
tune_grid(resamples = house_cv_folds,
grid = grid_knn,
control = control_grid(verbose = TRUE, allow_par = TRUE,
pkgs = c("outliers", "tidyverse")))
# Finalizamos clusters
stopCluster(make_cluster)
registerDoSEQ()
house_knn_vf |> collect_metrics()
# A tibble: 120 × 9
k weight dist .metric .estimator mean n std_err .config
<dbl> <chr> <dbl> <chr> <chr> <dbl> <int> <dbl> <chr>
1 8 gauss… 0.01 rmse standard 3.83e+4 32 1.19e+3 Prepro…
2 8 gauss… 0.01 rsq standard 7.71e-1 32 9.62e-3 Prepro…
3 20 gauss… 0.01 rmse standard 3.93e+4 32 1.29e+3 Prepro…
4 20 gauss… 0.01 rsq standard 7.64e-1 32 1.03e-2 Prepro…
5 60 gauss… 0.01 rmse standard 4.23e+4 32 1.36e+3 Prepro…
6 60 gauss… 0.01 rsq standard 7.39e-1 32 1.03e-2 Prepro…
7 80 gauss… 0.01 rmse standard 4.34e+4 32 1.37e+3 Prepro…
8 80 gauss… 0.01 rsq standard 7.30e-1 32 9.94e-3 Prepro…
9 100 gauss… 0.01 rmse standard 4.45e+4 32 1.39e+3 Prepro…
10 100 gauss… 0.01 rsq standard 7.19e-1 32 9.94e-3 Prepro…
# … with 110 more rows
Con show_best() se muestran a continuación los mejores modelos según `rsq``.
# Se muestran los mejores según rsq
house_knn_vf |> show_best("rsq")
# A tibble: 5 × 9
k weight dist .metric .estimator mean n std_err .config
<dbl> <chr> <dbl> <chr> <chr> <dbl> <int> <dbl> <chr>
1 8 gaussian 0.9 rsq standard 0.823 32 0.00774 Preproc…
2 8 inv 0.9 rsq standard 0.820 32 0.00802 Preproc…
3 20 gaussian 0.9 rsq standard 0.817 32 0.00825 Preproc…
4 20 inv 0.9 rsq standard 0.810 32 0.00872 Preproc…
5 60 gaussian 0.9 rsq standard 0.803 32 0.00825 Preproc…
Tras elegir el mejor, finalizamos el flujo con ese modelo elegido empleando la función finalize_workflow().
# Finalizamos flujo con el mejor modelo según rsq
best_house_knn_vf <-
house_knn_vf |> select_best("rsq")
final_wflow_knn <-
house_wflow_knn |>
finalize_workflow(best_house_knn_vf)
final_wflow_knn
══ Workflow ══════════════════════════════════════════════════════════
Preprocessor: Recipe
Model: nearest_neighbor()
── Preprocessor ──────────────────────────────────────────────────────
28 Recipe Steps
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• ...
• and 18 more steps.
── Model ─────────────────────────────────────────────────────────────
K-Nearest Neighbor Model Specification (regression)
Main Arguments:
neighbors = 8
weight_func = gaussian
dist_power = 0.9
Computational engine: kknn
A continuación, procederemos con el modelo para la receta de árboles.
El modelo consistirá en un árbol en modo regresión que implementaremos con decision_tree() y mode = "regression".
De nuevo, crearemos el modelo con tune() para poder definir posteriormente los parámetros de forma manual.
Generamos un flujo de trabajo en función del criterios de impureza de Gini.
Por otro lado, añadiremos también directamente el parámetro cost_complexity, que se usa en la fase de poda (prune), y que introduce un componente penalizador para evitar modelos sobreajustados.
# Modelo con tres parámetros tuneados
decision_tree_gini <-
decision_tree(tree_depth = tune("depth"), min_n = tune("min_n"),
cost_complexity = tune("cost")) |>
set_engine("rpart") |>
set_mode("regression")
# Flujo de trabajo
house_tree_gini_wflow <-
workflow() |>
add_recipe(house_rec_arbol) |>
add_model(decision_tree_gini)
Para este modelo, volvemos a evaluar a través del método de Validación Cruzada V-Folds, que nos proporcionará además una métrica media con la desviación típica para cada modelo.
Para ello, en primer lugar crearemos un grid con la función expand_grid que nos permitirán definir manualmente los parámetros que queramos sin estar condicionados por una función automática.
De nuevo, tras comprobar varias combinaciones, nos decantaremos finalmente por el siguiente.
grid_tree_gini <-
expand_grid("depth" = c(4, 6, 7, 8, 10),
"min_n" = c(45, 50, 55, 100, 150),
"cost" = 10^c(-5, -4, -3, -0.5))
Crearemos en este caso 100 modelos (\(5*5*4\)).
Aplicamos en esta fase el proceso de Validación Cruzada V-Folds.
set.seed(05492)
clusters <- detectCores() - 1
make_cluster <- makeCluster(clusters)
registerDoParallel(make_cluster)
house_tree_gini <-
house_tree_gini_wflow |>
tune_grid(resamples = house_cv_folds,
grid = grid_tree_gini,
control = control_grid(verbose = TRUE, allow_par = TRUE,
pkgs = c("outliers", "tidyverse"), save_pred = TRUE))
# Finalizamos clusters
stopCluster(make_cluster)
registerDoSEQ()
house_tree_gini |> collect_metrics()
# A tibble: 200 × 9
cost depth min_n .metric .estimator mean n std_err
<dbl> <dbl> <dbl> <chr> <chr> <dbl> <int> <dbl>
1 0.00001 4 45 rmse standard 47332. 32 1049.
2 0.00001 4 45 rsq standard 0.645 32 0.0110
3 0.0001 4 45 rmse standard 47332. 32 1049.
4 0.0001 4 45 rsq standard 0.645 32 0.0110
5 0.001 4 45 rmse standard 47331. 32 1049.
6 0.001 4 45 rsq standard 0.645 32 0.0110
7 0.316 4 45 rmse standard 62793. 32 1145.
8 0.316 4 45 rsq standard 0.372 32 0.0140
9 0.00001 4 50 rmse standard 46869. 32 1058.
10 0.00001 4 50 rsq standard 0.651 32 0.0107
# … with 190 more rows, and 1 more variable: .config <chr>
Como sucedía en el modelo KNN, procedemos a ver los mejores con show_best().
# Se muestran los mejores según rsq
house_tree_gini |> show_best("rsq", n = 10)
# A tibble: 10 × 9
cost depth min_n .metric .estimator mean n std_err .config
<dbl> <dbl> <dbl> <chr> <chr> <dbl> <int> <dbl> <chr>
1 0.00001 8 50 rsq standard 0.699 32 0.0105 Preproc…
2 0.0001 8 50 rsq standard 0.699 32 0.0105 Preproc…
3 0.001 8 50 rsq standard 0.699 32 0.0106 Preproc…
4 0.001 10 50 rsq standard 0.698 32 0.0107 Preproc…
5 0.00001 10 50 rsq standard 0.698 32 0.0107 Preproc…
6 0.0001 10 50 rsq standard 0.698 32 0.0107 Preproc…
7 0.00001 8 55 rsq standard 0.697 32 0.0104 Preproc…
8 0.0001 8 55 rsq standard 0.697 32 0.0104 Preproc…
9 0.0001 10 55 rsq standard 0.696 32 0.0106 Preproc…
10 0.00001 10 55 rsq standard 0.696 32 0.0106 Preproc…
Tras elegir el mejor, finalizamos el flujo con ese modelo elegido empleando la función finalize_workflow().
# Finalizamos flujo con el mejor modelo según rsq
best_house_tree_gini <-
house_tree_gini |> select_best("rsq")
final_wf_tree_gini <-
house_tree_gini_wflow |>
finalize_workflow(best_house_tree_gini)
final_wf_tree_gini
══ Workflow ══════════════════════════════════════════════════════════
Preprocessor: Recipe
Model: decision_tree()
── Preprocessor ──────────────────────────────────────────────────────
24 Recipe Steps
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• ...
• and 14 more steps.
── Model ─────────────────────────────────────────────────────────────
Decision Tree Model Specification (regression)
Main Arguments:
cost_complexity = 1e-05
tree_depth = 8
min_n = 50
Computational engine: rpart
A continuación un resumen de sus métricas para hacer la comparativa:
# Modelo KNN
house_knn_vf |>
collect_metrics() |>
filter(.config == pull(select(best_house_knn_vf, .config)))
# A tibble: 2 × 9
k weight dist .metric .estimator mean n std_err .config
<dbl> <chr> <dbl> <chr> <chr> <dbl> <int> <dbl> <chr>
1 8 gaussi… 0.9 rmse standard 3.40e+4 32 1.06e+3 Prepro…
2 8 gaussi… 0.9 rsq standard 8.23e-1 32 7.74e-3 Prepro…
# Modelo Árbol Gini
house_tree_gini |>
collect_metrics() |>
filter(.config == pull(select(best_house_tree_gini, .config)))
# A tibble: 2 × 9
cost depth min_n .metric .estimator mean n std_err .config
<dbl> <dbl> <dbl> <chr> <chr> <dbl> <int> <dbl> <chr>
1 0.00001 8 50 rmse standard 4.36e+4 32 1.16e+3 Prepro…
2 0.00001 8 50 rsq standard 6.99e-1 32 1.05e-2 Prepro…
Tanto en RMSE como en RSQ el modelo KNN parece mejor: presenta un menor RMSE y también un mayor RSQ.
Nos quedamos finalmente con el KNN.
# Ajuste
house_knn_fit <-
fit(final_wflow_knn, house_train)
# Predicciones
pred_knn <-
predict(house_knn_fit, house_test)
summary(pred_knn)
# Visualización de las predicciones para cada `Id`
final_pred <-
data.frame(Id = c(1461:2919), pred_knn) |>
dplyr::rename(SalePrice = .pred)
head(final_pred)
# Exportamos nuestro dataframe para subirlo a Kaggle
write.csv(final_pred, file = 'house_entrega_knn.csv', row.names = FALSE)
El tercer modelo que vamos a probar para predecir será la regresión univariante. En este tipo de regresión, seleccionaremos la variable que mayor correlación presenta respecto de nuestra variable objetivo SalePrice. A continuación, visualizaremos de nuevo la matriz de correlaciones y seleccionaremos la susodicha variable.
cor_matrix |>
corrplot(method = "number", tl.cex = 0.55, number.cex = 0.7, type = "lower")
A la vista de la matriz, la variable que mayor correlación presenta respecto de nuestra objetivo es GrLivArea con un coeficiente de correlación del 0.71. De esta manera, pasaremos a predecir nuestra objetivo SalePrice en función de GrLivArea.
Como la regresión univariante no admite ausentes ni en la predictora ni en la objetivo, lo que vamos a hacer es utilizar únicamente el conjunto de train en el que nuestra variable objetivo está completamente etiquetada.
ausentes <-
apply(house_train, 2, function(x) sum(is.na(x)))
ausentes_tb <-
tibble(Variable = names(house_train), Ausentes = ausentes) |>
filter(Variable == "GrLivArea" | Variable == "SalePrice")
ausentes_tb
# A tibble: 2 × 2
Variable Ausentes
<chr> <int>
1 GrLivArea 0
2 SalePrice 0
Como se puede observar, no hay ausentes ni en nuestra predictora ni en nuestra objetivo.
house[!is.na(house$SalePrice), ] |>
ggplot(aes(x = GrLivArea, y = SalePrice)) +
geom_point(col = "#EB9891") +
geom_smooth(method = "lm", se = FALSE, color = "black", aes(group = 1)) +
scale_y_log10(breaks= seq(0, 800000, by = 100000),
labels = comma_format(big.mark = " ", decimal.mark = ",")) +
scale_x_log10(labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_text_repel(aes(label = ifelse(house$GrLivArea[!is.na(house$SalePrice)] > 4000,
rownames(house), ""))) +
geom_text_repel(aes(label = ifelse(house$SalePrice[!is.na(house$SalePrice)] < 50000,
rownames(house), ""))) +
theme_minimal()
Si se grafica la relación siguiendo una escala logarítmica, podremos visualizar claramente esa correlación directa entre variables: a mayor superficie habitable, mayor precio. Además, podemos visualizar también varios valores atípicos. Respecto de la variableGrLivArea, tenemos las tuplas correspondientes a los Id 692, 1183, 524 y 1299. Pasaremos ahora a detectarlos numéricamente, pero antes transformaremos nuestras variables a logarítmicas.
house_log <-
house_train |>
mutate(SalePrice = log(SalePrice),
GrLivArea = log(GrLivArea))
cor_matrix <-
house_log |>
select(where(is.numeric)) |> cor(use = "pairwise.complete.obs", method = "pearson")
cor_matrix |>
corrplot(method = "number", tl.cex = 0.55, number.cex = 0.7, type = "lower")
El coeficiente de correlación ha aumentado ligeramente, de 0.71 a 0.73.
gg1 <- ggplot(house_train, aes(sample = SalePrice))+
stat_qq()+
stat_qq_line()+
labs(x = "SalePrice (lineal)", y = "") +
theme_minimal()
gg2 <- ggplot(house_train, aes(sample = GrLivArea))+
stat_qq()+
stat_qq_line()+
labs(x = "GrLivArea (lineal)", y = "") +
theme_minimal()
gg3 <- ggplot(house_log, aes(sample = SalePrice))+
stat_qq()+
stat_qq_line()+
labs(x = "SalePrice (logarítmica)", y = "") +
theme_minimal()
gg4 <- ggplot(house_log, aes(sample = GrLivArea))+
stat_qq()+
stat_qq_line()+
labs(x = "GrLivArea (logarítmica)", y = "") +
theme_minimal()
multiplot(gg1, gg2, gg3, gg4, cols = 2)
En ambas variables podemos observar como la transformación a logaritmo favorece la relación lineal. Ahora sí, pasemos a detectar outliers sobre nuestra predictora:
house_log |>
ggplot(aes(x = SalePrice)) +
geom_density(alpha = .8, fill="#EB9891") +
labs(title = "Distribución del precio de las viviendas (escala logarítmica)", x = "Precio", y = NULL) +
theme_minimal() +
theme(text = element_text(face = "bold", size = 15), plot.title = element_text(hjust = 0.5),
axis.title.x = element_text(vjust=-0.5)) +
scale_x_continuous(labels = comma_format(big.mark = " ", decimal.mark = ",")) +
scale_y_continuous(labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_vline(aes(xintercept = mean(SalePrice, na.rm = T), linetype = "Media"), colour = "black", size = .8) +
geom_vline(aes(xintercept = median(SalePrice, na.rm = T), linetype = "Mediana"), colour = "black", size = .8) +
scale_linetype_manual(name = "Medidas", values = c(Media = "solid", Mediana = "dotted"))
house_log |>
ggplot(aes(x = GrLivArea)) +
geom_density(alpha = .8, fill="#EB9891") +
labs(title = "Distribución de la superficie habitable de las viviendas en pies cuadrados (escala logarítmica)", x = "Pies cuadrados", y = NULL) +
theme_minimal() +
theme(text = element_text(face = "bold", size = 15), plot.title = element_text(hjust = 0.5),
axis.title.x = element_text(vjust=-0.5)) +
scale_x_continuous(labels = comma_format(big.mark = " ", decimal.mark = ",")) +
scale_y_continuous(labels = comma_format(big.mark = " ", decimal.mark = ",")) +
geom_vline(aes(xintercept = mean(GrLivArea, na.rm = T), linetype = "Media"), colour = "black", size = .8) +
geom_vline(aes(xintercept = median(GrLivArea, na.rm = T), linetype = "Mediana"), colour = "black", size = .8) +
scale_linetype_manual(name = "Medidas", values = c(Media = "solid", Mediana = "dotted"))
Ambas variables parecen seguir distribuciones más o menos simétricas, por lo que detectaremos los outliers de la predictora por la media.
# A tibble: 15 × 40
MSSubClass MSZoning LotFrontage LotArea Alley LandContour LandSlope
<fct> <fct> <dbl> <dbl> <fct> <fct> <fct>
1 30 RM 60 6324 None Lvl Gtl
2 75 RM 90 22950 None Lvl Gtl
3 75 RM 87 18386 None Lvl Gtl
4 60 RL 130 40094 None Bnk Gtl
5 30 RL 58 9098 None Lvl Gtl
6 20 RL 50 5000 None Low Mod
7 190 RH 60 10896 Pave Bnk Gtl
8 60 RL 104 21535 None Lvl Gtl
9 30 RM 50 6000 None Lvl Gtl
10 20 C (all) 50 9000 None Lvl Gtl
11 30 RL 60 8400 None Bnk Gtl
12 60 RL 118 35760 None Lvl Gtl
13 60 RL 160 15623 None Lvl Gtl
14 50 RL NA 14100 None Lvl Mod
15 60 RL 313 63887 None Bnk Gtl
# … with 33 more variables: Neighborhood <fct>, BldgType <fct>,
# HouseStyle <fct>, OverallCond <dbl>, YearBuilt <dbl>,
# RoofStyle <fct>, Exterior1st <fct>, ExterCond <int>,
# BsmtQual <int>, BsmtCond <int>, HeatingQC <int>,
# CentralAir <fct>, Electrical <fct>, `1stFlrSF` <dbl>,
# `2ndFlrSF` <dbl>, GrLivArea <dbl>, KitchenQual <dbl>,
# Functional <fct>, FireplaceQu <int>, GarageCars <dbl>, …
Durante la creación de la receta le imputaremos a estos outliers la media.
house_rec_uni <-
recipe(data = house_train, SalePrice ~ GrLivArea) |>
step_log(SalePrice, base = 10) |>
step_log(GrLivArea, base = 10) |>
step_mutate(GrLivArea = ifelse(abs(scores(GrLivArea, type = "z")) > 2.5, NA, GrLivArea)) |>
step_impute_mean(GrLivArea)
house_rec_uni
Recipe
Inputs:
role #variables
outcome 1
predictor 1
Operations:
Log transformation on SalePrice
Log transformation on GrLivArea
Variable mutation for ifelse(abs(scores(GrLivArea, typ...
Mean imputation for GrLivArea
Por último, horneamos nuestra receta para comprobar que todo funciona correctamente.
# A tibble: 1,460 × 2
GrLivArea SalePrice
<dbl> <dbl>
1 3.23 5.32
2 3.10 5.26
3 3.25 5.35
4 3.23 5.15
5 3.34 5.40
6 3.13 5.16
7 3.23 5.49
8 3.32 5.30
9 3.25 5.11
10 3.03 5.07
# … with 1,450 more rows
# Creación del modelo
reg_lineal <-
linear_reg() |> set_mode("regression") |> set_engine("lm")
# Creación del flujo
reg_wflow_uni <-
workflow() |>
add_model(reg_lineal) |>
add_recipe(house_rec_uni)
reg_wflow_uni
══ Workflow ══════════════════════════════════════════════════════════
Preprocessor: Recipe
Model: linear_reg()
── Preprocessor ──────────────────────────────────────────────────────
4 Recipe Steps
• step_log()
• step_log()
• step_mutate()
• step_impute_mean()
── Model ─────────────────────────────────────────────────────────────
Linear Regression Model Specification (regression)
Computational engine: lm
set.seed(05492)
reg_fit_uni <-
reg_wflow_uni |> fit(data = house_train)
reg_fit_uni
══ Workflow [trained] ════════════════════════════════════════════════
Preprocessor: Recipe
Model: linear_reg()
── Preprocessor ──────────────────────────────────────────────────────
4 Recipe Steps
• step_log()
• step_log()
• step_mutate()
• step_impute_mean()
── Model ─────────────────────────────────────────────────────────────
Call:
stats::lm(formula = ..y ~ ., data = data)
Coefficients:
(Intercept) GrLivArea
2.4439 0.8804
# Resumen del ajuste
tidy(reg_fit_uni)
# A tibble: 2 × 5
term estimate std.error statistic p.value
<chr> <dbl> <dbl> <dbl> <dbl>
1 (Intercept) 2.44 0.0752 32.5 1.27e-174
2 GrLivArea 0.880 0.0238 37.0 1.42e-211
Como se puede observar en las métricas, el modelo que se ha decidido crear se trata de un modelo doble logarítmico. La interpretación de \(\beta_1\) es la siguiente: si el logaritmo de la superficie habitable de una vivienda aumenta en un 1 % se estima que, en media, el logaritmo del precio de la vivienda se incrementará en un 0,88 %. Esta relación no es más que la elasticidad entre la superficie habitable de una vivienda y su precio y se trata, además, de una relación constante. Por otro lado, \(\beta_0\) puede interpretarse como el valor medio del precio de la vivienda cuando la superficie habitable toma valor cero. En este caso, \(\beta_0\) toma el valor de 2.44, aunque en muchas situaciones su interpretación carece de sentido.
Por otro lado, en la salida del modelo podemos observar como el p-valor coligado al parámetro \(\beta_1\) y al estadístico para este contraste es p = 1.415987e-211. Si trabajamos con un nivel de significación habitual del 5 % (0.05), el p-valor es menor que ese nivel de significación, por lo que podemos rechazar la hipótesis nula, esto es, la no significación individual de la variable log(GrLivArea) al 5 % de significación. Por ende, a la luz de los resultados, \(\beta_1\) es significativo individualmente, o lo que es lo mismo: el efecto de la variable log(GrLivArea) depende de la variable log(SalePrice) al 5 % de nivel de significación.
# Intervalo de confianza al 95 %
confint(reg_fit_uni |> extract_fit_engine())
2.5 % 97.5 %
(Intercept) 2.2963505 2.591491
GrLivArea 0.8336623 0.927106
Los intervalos de confianza nos ofrecen un rango de valores para \(\beta_j\) en donde no se rechazaría la hipótesis \(H_0\). En este caso concreto, el intervalo de confianza nos proporcionará un rango de valores para la variable explicativa log(GrLivArea) en donde podríamos comprobar su significación individual dentro del modelo.
reg_fit_uni |> extract_fit_engine() |> performance()
# Indices of model performance
AIC | AICc | BIC | R2 | R2 (adj.) | RMSE | Sigma
---------------------------------------------------------------------
-1931.953 | -1931.937 | -1916.095 | 0.484 | 0.483 | 0.125 | 0.125
check_model(reg_fit_uni |> extract_fit_engine())
En el gráfico, la línea es bastante atendencial, se ve bastante monótona en todo su recorrido. Comprobémoslo analíticamente con un ANOVA entre residuos y predictora:
adjustment <-
reg_fit_uni |> extract_fit_engine()
lm(adjustment$residuals ~ adjustment$fitted.values) |> anova()
Analysis of Variance Table
Response: adjustment$residuals
Df Sum Sq Mean Sq F value Pr(>F)
adjustment$fitted.values 1 0.000 0.000000 0 1
Residuals 1458 22.668 0.015547
Analysis of Variance Table
Response: adjustment$residuals
Df Sum Sq Mean Sq F value Pr(>F)
I(adjustment$fitted.values^2) 1 0.0000 0.000007 0.0004 0.9831
adjustment$fitted.values 1 0.0332 0.033153 2.1340 0.1443
Residuals 1457 22.6347 0.015535
Como se puede observar, el p-valor de la prueba F de Fisher-Snedecor es igual a 1, por lo que, a un nivel de significación del 0.05, no podemos rechazar la hipótesis nula, esto es, la no existencia de una relación lineal entre la predictora (GrLivArea) y sus residuos. Por tanto, ello implica que no parece existir tendencia lineal entre residuos y predictora en nuestro modelo. Lo mismo ocurre en el segundo caso: tampoco parece existir tendencia cuadrática entre residuos y predictora.
En el gráfico, la línea es bastante atendencial, se ve bastante monótona en todo su recorrido (excepto quizá en los extremos). Comprobémoslo analíticamente con un test de heterocedasticidad:
check_heteroscedasticity(adjustment)
Warning: Heteroscedasticity (non-constant error variance) detected (p < .001).
Según el test de heterocedasticidad, el supuesto de homocedasticidad no se cumple para nuestro modelo. Esto, muy a menudo, resulta poco realista, pues el supuesto implica que la variabilidad del error de la variable objetivo es la misma para cualquier nivel de nuestra predictora.
ggplot(
tibble("Observaciones" = 1:length(adjustment$residuals),
"Residuos" = adjustment$residuals),
aes(x = Observaciones, y = Residuos)) +
geom_point() +
geom_smooth(method = "lm", se = FALSE) +
theme_minimal()
En el gráfico, la recta de regresión se aprecia constante, en torno a 0, y los residuos se agrupan en torno a una banda también relativamente constante.
En el gráfico, los percentiles empíricos de los residuos se acoplan bastante bien a la diagonal, excepto en los extremos. Comprobémoslo analíticamente con un contraste de normalidad:
library(olsrr)
ols_test_normality(adjustment$residuals)
-----------------------------------------------
Test Statistic pvalue
-----------------------------------------------
Shapiro-Wilk 0.9642 0.0000
Kolmogorov-Smirnov 0.0685 0.0000
Cramer-von Mises 383.4014 0.0000
Anderson-Darling 11.4357 0.0000
-----------------------------------------------
Como se puede observar, los p-valor de los distintos test son iguales a 0, por lo que, a un nivel de significación del 0.05, podemos rechazar la hipótesis de normalidad. Por tanto, no podemos asumir a normalidad en nuestro modelo. Se han hecho pruebas con cambios en la receta añadiendo distintas transformaciones a nuestra objetivo (step_YeoJohnson(), step_BoxCox(), step_sqrt()), pero tampoco ha sido posible.
library(car)
durbinWatsonTest(adjustment)
lag Autocorrelation D-W Statistic p-value
1 0.01831851 1.963326 0.488
Alternative hypothesis: rho != 0
Como se puede observar, el p-valor de la prueba Durbin-Watson es igual a 0.476, por lo que, a un nivel de significación del 0.05, no podemos rechazar la hipótesis nula, esto es, la incorrelación de los residuos.
En el gráfico original de la diagnosis, todas las observaciones parecen estar contenidas entre las líneas punteadas al 0.5.
glance(reg_fit_uni)
# A tibble: 1 × 12
r.squared adj.r.squared sigma statistic p.value df logLik
<dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
1 0.484 0.483 0.125 1366. 1.42e-211 1 969.
# … with 5 more variables: AIC <dbl>, BIC <dbl>, deviance <dbl>,
# df.residual <int>, nobs <int>
Nuestro \(R^2\) es igual a 0.4838, por lo que el ratio de información explicada del modelo es del 48.38 %. Al ser un modelo tan simple con una única variable, la interpretación del resto de parámetros no tiene mucho sentido. Los dejaremos para los próximos modelos. A continuación, predeciremos y evaluaremos en un nuevo conjunto test que extraeremos de house_train, puesto que nuestro conjunto house_test original tiene la variable objetivo SalePrice repleta de NA.
# Hacemos el split. Lo hacemos del 0.6 ya que hay pocos datos en train.
set.seed(05492)
split_house <-
initial_split(house_train, prop = 0.6)
# Predecimos en test
fit_reg_uni <-
reg_wflow_uni |> last_fit(split = split_house)
# Evaluamos en test
fit_reg_uni |> collect_metrics()
# A tibble: 2 × 4
.metric .estimator .estimate .config
<chr> <chr> <dbl> <chr>
1 rmse standard 0.126 Preprocessor1_Model1
2 rsq standard 0.455 Preprocessor1_Model1
# Errores en test
pred_test <-
fit_reg_uni |>
collect_predictions() |>
mutate(error = SalePrice - .pred)
pred_test
# A tibble: 584 × 6
id .pred .row SalePrice .config error
<chr> <dbl> <int> <dbl> <chr> <dbl>
1 train/test split 5.17 2 5.26 Preprocessor1_Model1 0.0874
2 train/test split 5.31 3 5.35 Preprocessor1_Model1 0.0394
3 train/test split 5.29 7 5.49 Preprocessor1_Model1 0.198
4 train/test split 5.37 8 5.30 Preprocessor1_Model1 -0.0715
5 train/test split 5.31 9 5.11 Preprocessor1_Model1 -0.194
6 train/test split 5.09 11 5.11 Preprocessor1_Model1 0.0179
7 train/test split 5.04 13 5.16 Preprocessor1_Model1 0.116
8 train/test split 5.17 15 5.20 Preprocessor1_Model1 0.0273
9 train/test split 5.08 17 5.17 Preprocessor1_Model1 0.0929
10 train/test split 5.20 20 5.14 Preprocessor1_Model1 -0.0520
# … with 574 more rows
En esta tabla se puede observar cómo, en nuestro modelo doble logarítmico, las estimaciones se ajustan en gran medida a los valores reales de la objetivo SalePrice. Visualizaremos ahora estos valores a través de un gráfico.
g1 <- pred_test |>
ggplot(mapping = aes(x = .pred, y = SalePrice)) +
geom_point(color = "#56BCC2", alpha = 0.6, size = 4) +
geom_abline(intercept = 0, slope = 1, color = "#EB9891", size = 1.2) +
theme_minimal() +
labs(title = "Resultados de la regresión lineal univariante",
subtitle = "Los valores predichos deberían estar cercanos a la diagonal",
x = "Predicciones",
y = "Valores reales")
g2 <- pred_test |>
select(.pred, SalePrice) |>
gather(Distribución, value) |>
ggplot(aes(x = value, color = Distribución, fill = Distribución)) +
geom_density(alpha = 0.6) +
theme_minimal() +
labs(title = "Distribución de las predicciones sobre los valores reales de SalePrice",
x = "Distribuciones",
y = "Frecuencia")
multiplot(g1, g2)
Como se puede apreciar en los gráficos, las predicciones son un poco reguleras, principalmente por el reducido tamaño muestral y por el bajo \(R^2\) (finalmente, el ratio de información explicada del modelo fue del 49.88 %). Tampoco se han cumplido varias de las hipótesis. Desde una regresión lineal univariante, en la que únicamente disponemos de una única variable para predecir la objetivo, poco más se puede hacer. En los siguientes apartados avanzaremos hacia la regresión multivariante con la creación de dos nuevos modelos.
Para este nuevo modelo de regresión vamos a emplear todas las variables. Volvemos a generar la receta, pero esta vez enfrentando nuestra objetivo al total de variables predictoras que tenemos en el dataset.
Volvemos a definir una nueva receta indicándole el conjunto donde tenemos validación y train, y enfrentando nuestra variable objetivo SalePrice a todas las demás.
Después, asignamos posibles roles, sujetos a modificación, que nos permitan diferenciar acciones entre las variables (sobre todo en la sección outliers).
# Receta
house_rec_multi <-
# Fórmula y datos
recipe(data = house_train, SalePrice ~ .)|>
# Roles
add_role(where(is.factor),
new_role = "cualitativa") |>
add_role(where(is.numeric),
new_role = "cuantitativa") |>
add_role(RemodAdd, New,
new_role = "binaria") |>
add_role(LotFrontage, LotArea, `1stFlrSF`, `2ndFlrSF`, GrLivArea, Age,
new_role = "mediana") |>
add_role(all_nominal_predictors(),
new_role = "moda") |>
add_role(ExterCond, BsmtCond, HeatingQC, GarageQual, KitchenQual,
new_role = "imputar mediana") |>
add_role(OverallCond, BsmtQual, FireplaceQu, TotBath, YearBuilt, GarageCars,
new_role = "imputar media") |>
add_role(SalePrice, GrLivArea,`1stFlrSF`, `2ndFlrSF`,
new_role = "log")
En este apartado, reagrupamos las variables cualitativas con step_mutate con lo definido en la fase de exploración. Vamos a reducir el número de categorías para que, al dummyficar, no se creen en exceso sin ser realmente necesarias.
house_rec_multi <-
house_rec_multi |>
step_mutate(Alley =
fct_collapse(Alley,
Grvl_Pave = c("Grvl", "Pave"))) |>
step_mutate(LandSlope =
fct_collapse(LandSlope,
Mod_Sev = c("Mod", "Sev"))) |>
step_mutate(LandContour =
fct_collapse(LandContour,
Bnk_HLS_Low = c("Bnk", "HLS", "Low"))) |>
step_mutate(RoofStyle =
fct_collapse(RoofStyle,
Hip = c("Flat", "Shed", "Mansard", "Hip"),
Gam_Gable = c("Gambrel", "Gable"))) |>
step_mutate(Exterior1st =
fct_collapse(Exterior1st,
`<150k` = c("AsphShn", "BrkComm", "BrkFace",
"CemntBd", "HdBoard", "MetalSd",
"Stucco", "Wd Sdng", "WdShing",
"AsbShng"),
`>150k` = c("CBlock", "ImStucc", "Plywood",
"Stone", "VinylSd"))) |>
step_mutate(PavedDrive =
fct_collapse(PavedDrive,
N_P = c("N", "P"))) |>
step_mutate(Fence =
fct_collapse(Fence,
MnPrv = c("GdWo", "MnWw"),
None = c("GdPrv"))) |>
step_mutate(MSSubClass =
fct_collapse(MSSubClass,
`<150k` = c("30", "40", "45", "50",
"80", "85", "120", "160",
"180", "190"),
`>150k` = c("20", "60", "70", "75",
"90", "150"))) |>
step_mutate(MSZoning =
fct_collapse(MSZoning,
RM = c("C (all)", "RH"),
RL = c("FV"))) |>
step_mutate(Neighborhood =
fct_collapse(Neighborhood,
`<125k` = c("BrDale", "BrkSide", "Edwards",
"IDOTRR", "MeadowV", "OldTown"),
`125k-150k` = c("Blmngtn", "NAmes", "NPkVill",
"SWISU", "Sawyer"),
`150k-200k` = c("Blueste", "CollgCr", "Gilbert",
"Mitchel", "NWAmes", "SawyerW"),
`>200k` = c("ClearCr", "Crawfor", "NoRidge",
"NridgHt", "Somerst", "StoneBr",
"Timber", "Veenker"))) |>
step_mutate(BldgType =
fct_collapse(BldgType,
Fam_TwnhsE = c("1Fam", "TwnhsE"),
fmCon_Dup_Twnhs = c("2fmCon", "Duplex", "Twnhs"))) |>
step_mutate(HouseStyle =
fct_collapse(HouseStyle,
Story_Fin_SLvl = c("2Story", "2.5Fin", "SLvl"),
fmCon_Dup_Twnhs = c("1.5Fin", "1.5Unf", "1Story",
"2.5Unf", "SFoyer"))) |>
step_mutate(Electrical =
fct_collapse(Electrical,
Other = c("FuseA", "FuseF", "FuseP", "Mix"))) |>
step_mutate(Functional =
fct_collapse(Functional,
Other = c("Maj1", "Maj2", "Min1", "Min2",
"Mod", "Sev"))) |>
step_mutate(SaleType =
fct_collapse(SaleType,
Con_CWD_New = c("Con", "CWD", "New"),
WD = c("COD", "ConLD", "ConLI", "ConLw",
"Oth"))) |>
step_mutate(SaleCondition =
fct_collapse(SaleCondition,
Normal_Others = c("Abnorml", "AdjLand", "Alloca",
"Family"))) |>
step_mutate(MoSold =
fct_collapse(MoSold,
Jan_Ap_May_Oct = c("1", "4", "5", "10"),
Mar_Jun_Jul = c("3", "6", "7"),
Feb_Aug_Sep_Nov_Dec = c("2", "8", "9",
"11", "12"))) |>
step_mutate(YrSold =
fct_collapse(YrSold,
`2006-2009` = c("2006", "2007", "2008", "2009")))
En los métodos de regresión, la recategorización de las variables cuantitativas no es estrictamente necesaria. Por tanto, en este receta, del total de variables cuantitativas no vamos a modificar ninguna.
multiplot(box1, box2, box3, box4, box5, box6, box7, box8, box9, box10, cols = 2)
multiplot(box11, box12, box13, box14, box15, box16, box17, cols = 2)
Si observamos estos gráficos de cajas y bigotes, todas nuestras variables cuantitativas continuas son asimétricas, por lo que se detectarán los outliers y se imputarán los ausentes por la mediana. Para el caso de las semi-cuantitativas (KitchenQual, por ejemplo) se imputarán directamente los ausentes por uno de los dos métodos en función de la simetría o asimetría de la variable. Para el resto de variables cualitativas, imputamos directamente por la moda.
house_rec_multi <-
house_rec_multi |>
# Detección de outliers por la mediana y por la media
step_mutate(across(has_role("mediana"), function(x) { ifelse(abs(scores(x, type = "mad")) > 3 & !is.na(x), NA, x) })) |>
# Imputación de ausentes por la mediana, la media y la moda
step_impute_median(has_role("mediana")) |>
step_impute_median(has_role("imputar_mediana")) |>
step_impute_mean(has_role("imputar_media")) |>
step_impute_mode(has_role("moda"))
A fin de que nuestro modelo cumpla los test de normalidad de los residuos, vamos a incluir en la receta la función step_YeoJohnson, que se emplea para eliminar la asimetría en los datos numéricos del modelo.
house_rec_multi <-
house_rec_multi |>
step_YeoJohnson(all_numeric_predictors(), -has_role("log"))
Transformamos algunas variables a logarítmicas para mejorar su relación.
Normalizamos nuestras variables por rango para que todas tengan el mismo peso, entre 0 y 1.
house_rec_multi <-
house_rec_multi |>
step_normalize(all_numeric_predictors())
Como esta receta es para un modelo de regresión multivariante, debemos dummyficar nuestras variables cualitativas. Para ello, tomamos todas las nominales, menos nuestra variable objetivo.
house_rec_multi <-
house_rec_multi |>
step_dummy(all_nominal_predictors())
house_rec_multi <-
house_rec_multi |>
step_zv(all_predictors())
Aplicamos el filtro de correlación a nuestras variables numéricas.
house_rec_multi <-
house_rec_multi |>
step_corr(all_numeric_predictors(), threshold = 0.6)
Por último, horneamos nuestra receta para comprobar que todas nuestras nuevas variables recategorizadas se hayan creado correctamente.
# A tibble: 1,460 × 35
LotArea OverallCond ExterCond BsmtQual BsmtCond HeatingQC
<dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
1 -0.141 -0.477 -0.194 0.576 -0.00502 0.930
2 0.106 2.01 -0.194 0.576 -0.00502 0.930
3 0.414 -0.477 -0.194 0.576 -0.00502 0.930
4 0.0955 -0.477 -0.194 -0.695 3.56 -0.330
5 0.877 -0.477 -0.194 0.576 -0.00502 0.930
6 0.857 -0.477 -0.194 0.576 -0.00502 0.930
7 0.201 -0.477 -0.194 2.13 -0.00502 0.930
8 0.257 0.440 -0.194 0.576 -0.00502 0.930
9 -0.761 -0.477 -0.194 -0.695 -0.00502 -0.330
10 -0.391 0.440 -0.194 -0.695 -0.00502 0.930
# … with 1,450 more rows, and 29 more variables: `1stFlrSF` <dbl>,
# `2ndFlrSF` <dbl>, GrLivArea <dbl>, FireplaceQu <dbl>,
# GarageCars <dbl>, GarageQual <dbl>, RemodAdd <dbl>, New <dbl>,
# Age <dbl>, SalePrice <dbl>, MSSubClass_X.150k <dbl>,
# MSZoning_RL <dbl>, Alley_None <dbl>, LandContour_Lvl <dbl>,
# LandSlope_Mod_Sev <dbl>, Neighborhood_X150k.200k <dbl>,
# Neighborhood_X.200k <dbl>, BldgType_fmCon_Dup_Twnhs <dbl>, …
# Creación del flujo
reg_wflow_multi <-
workflow() |>
add_model(reg_lineal) |>
add_recipe(house_rec_multi)
reg_wflow_multi
══ Workflow ══════════════════════════════════════════════════════════
Preprocessor: Recipe
Model: linear_reg()
── Preprocessor ──────────────────────────────────────────────────────
29 Recipe Steps
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• ...
• and 19 more steps.
── Model ─────────────────────────────────────────────────────────────
Linear Regression Model Specification (regression)
Computational engine: lm
set.seed(05492)
reg_fit_multi <-
reg_wflow_multi |> fit(data = house_train)
reg_fit_multi
══ Workflow [trained] ════════════════════════════════════════════════
Preprocessor: Recipe
Model: linear_reg()
── Preprocessor ──────────────────────────────────────────────────────
29 Recipe Steps
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• step_mutate()
• ...
• and 19 more steps.
── Model ─────────────────────────────────────────────────────────────
Call:
stats::lm(formula = ..y ~ ., data = data)
Coefficients:
(Intercept) LotArea
11.6420832 0.0186542
OverallCond ExterCond
0.0466561 -0.0011050
BsmtQual BsmtCond
0.0718343 0.0044886
HeatingQC `1stFlrSF`
0.0216881 0.0319568
`2ndFlrSF` GrLivArea
-0.0230072 0.1635733
FireplaceQu GarageCars
0.0260895 0.0492290
GarageQual RemodAdd
0.0161250 -0.0192380
New Age
0.0025791 -0.0418791
MSSubClass_X.150k MSZoning_RL
-0.0138932 0.0641965
Alley_None LandContour_Lvl
0.0292293 0.0263570
LandSlope_Mod_Sev Neighborhood_X150k.200k
0.0253077 0.0383263
Neighborhood_X.200k BldgType_fmCon_Dup_Twnhs
0.1641738 -0.0502749
RoofStyle_Gam_Gable Exterior1st_X.150k
-0.0199646 -0.0006763
CentralAir_Y Electrical_SBrkr
0.0688261 0.0015624
Functional_Typ PavedDrive_Y
0.1148388 0.0462485
Fence_MnPrv MoSold_Feb_Aug_Sep_Nov_Dec
-0.0090422 -0.0073745
MoSold_Mar_Jun_Jul YrSold_X2010
-0.0003882 0.0058141
SaleCondition_Normal
0.0460244
# Resumen del ajuste
tidy(reg_fit_multi)
# A tibble: 35 × 5
term estimate std.error statistic p.value
<chr> <dbl> <dbl> <dbl> <dbl>
1 (Intercept) 11.6 0.0351 332. 0
2 LotArea 0.0187 0.00521 3.58 3.50e- 4
3 OverallCond 0.0467 0.00469 9.94 1.49e-22
4 ExterCond -0.00110 0.00413 -0.268 7.89e- 1
5 BsmtQual 0.0718 0.00587 12.2 8.35e-33
6 BsmtCond 0.00449 0.00428 1.05 2.94e- 1
7 HeatingQC 0.0217 0.00479 4.52 6.59e- 6
8 `1stFlrSF` 0.0320 0.0140 2.28 2.27e- 2
9 `2ndFlrSF` -0.0230 0.0143 -1.61 1.08e- 1
10 GrLivArea 0.164 0.0166 9.87 2.95e-22
# … with 25 more rows
# A tibble: 6 × 5
term estimate std.error statistic p.value
<chr> <dbl> <dbl> <dbl> <dbl>
1 ExterCond -0.00110 0.00413 -0.268 0.789
2 New 0.00258 0.00459 0.561 0.575
3 Exterior1st_X.150k -0.000676 0.00915 -0.0739 0.941
4 Electrical_SBrkr 0.00156 0.0151 0.104 0.917
5 MoSold_Mar_Jun_Jul -0.000388 0.00882 -0.0440 0.965
6 YrSold_X2010 0.00581 0.0117 0.498 0.618
Como se puede observar en las métricas, los p-valor coligados a los parámetros \(\beta\) de las variables ExterCond, New, Exterior1st_X.150k, Electrical_SBrkr, MoSold_Mar_Jun_Jul y YrSold_X2010 son superiores a 0.05. Si trabajamos con un nivel de significación habitual del 5 % (0.05), el p-valor es mayor que ese nivel de significación, por lo que no podemos rechazar la hipótesis nula, esto es, la no significación individual de estas variables al 5 % de significación. Por ende, a la luz de los resultados, ExterCond, New, Exterior1st_X.150k, Electrical_SBrkr, MoSold_Mar_Jun_Jul y YrSold_X2010 son variables no significativas individualmente. Optaremos por eliminarlas posteriormente en la selección de modelos porque solo aportan ruido.
reg_fit_multi |> extract_fit_engine() |> performance()
# Indices of model performance
AIC | AICc | BIC | R2 | R2 (adj.) | RMSE | Sigma
---------------------------------------------------------------------
-1543.682 | -1541.810 | -1353.379 | 0.879 | 0.876 | 0.139 | 0.141
Nuestro \(R^2\) es igual a 0.879, por lo que el ratio de información explicada del modelo es del 88 %. Este modelo ha duplicado el ratio de información explicada respecto del modelo univariante anterior.
check_model(reg_fit_multi |> extract_fit_engine())
En el gráfico, la línea es bastante atendencial, se ve bastante monótona en todo su recorrido. Comprobémoslo analíticamente con un ANOVA entre residuos y predictora:
adjustment <-
reg_fit_multi |> extract_fit_engine()
lm(adjustment$residuals ~ adjustment$fitted.values) |> anova()
Analysis of Variance Table
Response: adjustment$residuals
Df Sum Sq Mean Sq F value Pr(>F)
adjustment$fitted.values 1 0.000 0.000000 0 1
Residuals 1458 28.267 0.019387
Analysis of Variance Table
Response: adjustment$residuals
Df Sum Sq Mean Sq F value Pr(>F)
I(adjustment$fitted.values^2) 1 0.0001 0.000118 0.0061 0.9376883
adjustment$fitted.values 1 0.2483 0.248296 12.9119 0.0003374
Residuals 1457 28.0183 0.019230
I(adjustment$fitted.values^2)
adjustment$fitted.values ***
Residuals
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Como se puede observar, el p-valor de la prueba F de Fisher-Snedecor es igual a 1, por lo que, a un nivel de significación del 0.05, no podemos rechazar la hipótesis nula, esto es, la no existencia de una relación lineal entre las predictoras y sus residuos. Por tanto, ello implica que no parece existir tendencia lineal entre residuos y predictoras en nuestro modelo. Lo mismo ocurre en el segundo caso: tampoco parece existir tendencia cuadrática entre residuos y predictoras.
En el gráfico, la línea sigue una fase decreciente al principio y creciente al final, se ve poco monótona en todo su recorrido. Comprobémoslo analíticamente con un test de heterocedasticidad:
check_heteroscedasticity(adjustment)
Warning: Heteroscedasticity (non-constant error variance) detected (p < .001).
Según el test de heterocedasticidad, el supuesto de homocedasticidad no se cumple para nuestro modelo. Esto, muy a menudo, resulta poco realista, pues el supuesto implica que la variabilidad del error de la variable objetivo es la misma para cualquier nivel de nuestra predictora.
ggplot(
tibble("Observaciones" = 1:length(adjustment$residuals),
"Residuos" = adjustment$residuals),
aes(x = Observaciones, y = Residuos)) +
geom_point() +
geom_smooth(method = "lm", se = FALSE) +
theme_minimal()
En el gráfico, la recta de regresión se aprecia constante, en torno a 0, y los residuos se agrupan en torno a una banda también relativamente constante.
En el gráfico original de la diagnosis, prácticamente todas las observaciones están en la franja verde. Tan solo tres están en la roja (es el mejor resultado que obtuvimos).
En el gráfico, los percentiles empíricos de los residuos se acoplan bastante bien a la diagonal, excepto en los extremos. Comprobémoslo analíticamente con un contraste de normalidad:
library(olsrr)
ols_test_normality(adjustment$residuals)
-----------------------------------------------
Test Statistic pvalue
-----------------------------------------------
Shapiro-Wilk 0.9313 0.0000
Kolmogorov-Smirnov 0.0659 0.0000
Cramer-von Mises 374.2219 0.0000
Anderson-Darling 11.2063 0.0000
-----------------------------------------------
Como se puede observar, los p-valor de los distintos test son iguales a 0, por lo que, a un nivel de significación del 0.05, podemos rechazar la hipótesis de normalidad. Por tanto, no podemos asumir a normalidad en nuestro modelo. Se han hecho pruebas con cambios en la receta añadiendo distintas transformaciones a nuestra objetivo (step_YeoJohnson(), step_BoxCox(), step_sqrt()), pero tampoco ha sido posible.
library(car)
durbinWatsonTest(adjustment)
lag Autocorrelation D-W Statistic p-value
1 0.007249536 1.985144 0.828
Alternative hypothesis: rho != 0
Como se puede observar, el p-valor de la prueba Durbin-Watson es igual a 0.772, por lo que, a un nivel de significación del 0.05, no podemos rechazar la hipótesis nula, esto es, la incorrelación de los residuos.
En el gráfico original de la diagnosis, todas las observaciones parecen estar contenidas entre las líneas punteadas al 0.5.
glance(reg_fit_multi)
# A tibble: 1 × 12
r.squared adj.r.squared sigma statistic p.value df logLik AIC
<dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
1 0.879 0.876 0.141 303. 0 34 808. -1544.
# … with 4 more variables: BIC <dbl>, deviance <dbl>,
# df.residual <int>, nobs <int>
Nuestro \(R^2\) es igual a 0.8786, por lo que el ratio de información explicada del modelo es del 88 %. Este modelo ha duplicado el ratio de información explicada respecto del modelo univariante anterior.
# Hacemos el split. Lo hacemos del 0.6 ya que hay pocos datos en train.
set.seed(05492)
split_house <-
initial_split(house_train, prop = 0.6)
# Predecimos en test
reg_fit_multi <-
reg_wflow_multi |> last_fit(split = split_house)
# Evaluamos en test
reg_fit_multi |> collect_metrics()
# A tibble: 2 × 4
.metric .estimator .estimate .config
<chr> <chr> <dbl> <chr>
1 rmse standard 0.152 Preprocessor1_Model1
2 rsq standard 0.850 Preprocessor1_Model1
set.seed(05492)
# Errores en test
pred_test <-
reg_fit_multi |>
collect_predictions() |>
mutate(error = SalePrice - .pred)
pred_test
# A tibble: 584 × 6
id .pred .row SalePrice .config error
<chr> <dbl> <int> <dbl> <chr> <dbl>
1 train/test split 12.3 2 12.1 Preprocessor1_Mode… -0.232
2 train/test split 12.2 3 12.3 Preprocessor1_Mode… 0.112
3 train/test split 12.6 7 12.6 Preprocessor1_Mode… 0.0661
4 train/test split 12.3 8 12.2 Preprocessor1_Mode… -0.132
5 train/test split 11.6 9 11.8 Preprocessor1_Mode… 0.137
6 train/test split 11.7 11 11.8 Preprocessor1_Mode… 0.0424
7 train/test split 11.6 13 11.9 Preprocessor1_Mode… 0.233
8 train/test split 11.8 15 12.0 Preprocessor1_Mode… 0.146
9 train/test split 11.9 17 11.9 Preprocessor1_Mode… -0.00848
10 train/test split 11.6 20 11.8 Preprocessor1_Mode… 0.195
# … with 574 more rows
En esta tabla se puede observar cómo, en nuestro modelo, las estimaciones se ajustan en gran medida a los valores reales de la objetivo SalePrice. Visualizaremos ahora estos valores a través de un gráfico.
g1 <- pred_test |>
ggplot(mapping = aes(x = .pred, y = SalePrice)) +
geom_point(color = "#56BCC2", alpha = 0.6, size = 4) +
geom_abline(intercept = 0, slope = 1, color = "#EB9891", size = 1.2) +
theme_minimal() +
labs(title = "Resultados de la regresión lineal multivariante (modelo saturado)",
subtitle = "Los valores predichos deberían estar cercanos a la diagonal",
x = "Predicciones",
y = "Valores reales")
g2 <- pred_test |>
select(.pred, SalePrice) |>
gather(Distribución, value) |>
ggplot(aes(x = value, color = Distribución, fill = Distribución)) +
geom_density(alpha = 0.6) +
theme_minimal() +
labs(title = "Distribución de las predicciones sobre los valores reales de SalePrice",
x = "Distribuciones",
y = "Frecuencia")
multiplot(g1, g2)
Como se puede apreciar en los gráficos, las predicciones son bastante buenas, bastante mejores que en el modelo univariante. A pesar de que se han incumplido algunos de los supuestos de regresión, la recta se ajusta bastante bien a los datos. A continuación, crearemos un modelo a caballo entre el univariante y el saturado, incluyendo únicamente en el modelo una selección del total de variables del dataset.
Para este último modelo de regresión vamos a emplear únicamente un conjunto de variables seleccionadas a través de distintos criterios de información.
Para la selección de modelos, pasaremos directamente la receta a lm() para proceder a la regresión contra las variables escogidas.
En este caso, y tras varias pruebas, hemos decidido emplear el criterio de información BIC, por lo que le aplicaremos a la función stepAIC su penalización correspondiente:
Start: AIC=-5503.97
SalePrice ~ LotArea + OverallCond + ExterCond + BsmtQual + BsmtCond +
HeatingQC + `1stFlrSF` + `2ndFlrSF` + GrLivArea + FireplaceQu +
GarageCars + GarageQual + RemodAdd + New + Age + MSSubClass_X.150k +
MSZoning_RL + Alley_None + LandContour_Lvl + LandSlope_Mod_Sev +
Neighborhood_X150k.200k + Neighborhood_X.200k + BldgType_fmCon_Dup_Twnhs +
RoofStyle_Gam_Gable + Exterior1st_X.150k + CentralAir_Y +
Electrical_SBrkr + Functional_Typ + PavedDrive_Y + Fence_MnPrv +
MoSold_Feb_Aug_Sep_Nov_Dec + MoSold_Mar_Jun_Jul + YrSold_X2010 +
SaleCondition_Normal
Df Sum of Sq RSS AIC
- MoSold_Mar_Jun_Jul 1 0.00004 28.267 -5511.2
- Exterior1st_X.150k 1 0.00011 28.267 -5511.2
- Electrical_SBrkr 1 0.00021 28.267 -5511.2
- ExterCond 1 0.00142 28.268 -5511.2
- YrSold_X2010 1 0.00493 28.272 -5511.0
- New 1 0.00625 28.273 -5510.9
- MoSold_Feb_Aug_Sep_Nov_Dec 1 0.01073 28.277 -5510.7
- Fence_MnPrv 1 0.01344 28.280 -5510.6
- BsmtCond 1 0.02185 28.288 -5510.1
- LandSlope_Mod_Sev 1 0.03171 28.298 -5509.6
- MSSubClass_X.150k 1 0.03942 28.306 -5509.2
- `2ndFlrSF` 1 0.05127 28.318 -5508.6
- Alley_None 1 0.06139 28.328 -5508.1
- LandContour_Lvl 1 0.06187 28.328 -5508.1
- RoofStyle_Gam_Gable 1 0.08215 28.349 -5507.0
- `1stFlrSF` 1 0.10321 28.370 -5505.9
<none> 28.267 -5504.0
- PavedDrive_Y 1 0.17295 28.440 -5502.3
- Neighborhood_X150k.200k 1 0.22259 28.489 -5499.8
- BldgType_fmCon_Dup_Twnhs 1 0.22948 28.496 -5499.4
- LotArea 1 0.25478 28.521 -5498.2
- GarageQual 1 0.27613 28.543 -5497.1
- CentralAir_Y 1 0.28574 28.552 -5496.6
- RemodAdd 1 0.33524 28.602 -5494.0
- SaleCondition_Normal 1 0.33567 28.602 -5494.0
- HeatingQC 1 0.40587 28.672 -5490.4
- MSZoning_RL 1 0.52610 28.793 -5484.3
- FireplaceQu 1 0.63386 28.901 -5478.9
- Age 1 0.93476 29.201 -5463.8
- Functional_Typ 1 1.04377 29.310 -5458.3
- GarageCars 1 1.57238 29.839 -5432.2
- GrLivArea 1 1.93106 30.198 -5414.8
- OverallCond 1 1.95981 30.227 -5413.4
- BsmtQual 1 2.96959 31.236 -5365.4
- Neighborhood_X.200k 1 3.12340 31.390 -5358.2
Step: AIC=-5511.25
SalePrice ~ LotArea + OverallCond + ExterCond + BsmtQual + BsmtCond +
HeatingQC + `1stFlrSF` + `2ndFlrSF` + GrLivArea + FireplaceQu +
GarageCars + GarageQual + RemodAdd + New + Age + MSSubClass_X.150k +
MSZoning_RL + Alley_None + LandContour_Lvl + LandSlope_Mod_Sev +
Neighborhood_X150k.200k + Neighborhood_X.200k + BldgType_fmCon_Dup_Twnhs +
RoofStyle_Gam_Gable + Exterior1st_X.150k + CentralAir_Y +
Electrical_SBrkr + Functional_Typ + PavedDrive_Y + Fence_MnPrv +
MoSold_Feb_Aug_Sep_Nov_Dec + YrSold_X2010 + SaleCondition_Normal
Df Sum of Sq RSS AIC
- Exterior1st_X.150k 1 0.00011 28.267 -5518.5
- Electrical_SBrkr 1 0.00021 28.267 -5518.5
- ExterCond 1 0.00141 28.268 -5518.5
- YrSold_X2010 1 0.00511 28.272 -5518.3
- New 1 0.00627 28.273 -5518.2
- Fence_MnPrv 1 0.01342 28.280 -5517.8
- MoSold_Feb_Aug_Sep_Nov_Dec 1 0.01343 28.280 -5517.8
- BsmtCond 1 0.02188 28.289 -5517.4
- LandSlope_Mod_Sev 1 0.03197 28.299 -5516.9
- MSSubClass_X.150k 1 0.03938 28.306 -5516.5
- `2ndFlrSF` 1 0.05123 28.318 -5515.9
- Alley_None 1 0.06164 28.328 -5515.4
- LandContour_Lvl 1 0.06209 28.329 -5515.3
- RoofStyle_Gam_Gable 1 0.08215 28.349 -5514.3
- `1stFlrSF` 1 0.10353 28.370 -5513.2
<none> 28.267 -5511.2
- PavedDrive_Y 1 0.17292 28.440 -5509.6
- Neighborhood_X150k.200k 1 0.22263 28.489 -5507.1
- BldgType_fmCon_Dup_Twnhs 1 0.22944 28.496 -5506.7
- LotArea 1 0.25492 28.522 -5505.4
- GarageQual 1 0.27609 28.543 -5504.3
- CentralAir_Y 1 0.28666 28.553 -5503.8
- RemodAdd 1 0.33528 28.602 -5501.3
- SaleCondition_Normal 1 0.33564 28.602 -5501.3
- HeatingQC 1 0.40645 28.673 -5497.7
- MSZoning_RL 1 0.52662 28.793 -5491.6
- FireplaceQu 1 0.63444 28.901 -5486.1
- Age 1 0.93473 29.201 -5471.0
- Functional_Typ 1 1.04385 29.311 -5465.6
- GarageCars 1 1.57360 29.840 -5439.4
- GrLivArea 1 1.93324 30.200 -5421.9
- OverallCond 1 1.96045 30.227 -5420.6
- BsmtQual 1 2.98093 31.248 -5372.2
- Neighborhood_X.200k 1 3.15004 31.417 -5364.3
Step: AIC=-5518.53
SalePrice ~ LotArea + OverallCond + ExterCond + BsmtQual + BsmtCond +
HeatingQC + `1stFlrSF` + `2ndFlrSF` + GrLivArea + FireplaceQu +
GarageCars + GarageQual + RemodAdd + New + Age + MSSubClass_X.150k +
MSZoning_RL + Alley_None + LandContour_Lvl + LandSlope_Mod_Sev +
Neighborhood_X150k.200k + Neighborhood_X.200k + BldgType_fmCon_Dup_Twnhs +
RoofStyle_Gam_Gable + CentralAir_Y + Electrical_SBrkr + Functional_Typ +
PavedDrive_Y + Fence_MnPrv + MoSold_Feb_Aug_Sep_Nov_Dec +
YrSold_X2010 + SaleCondition_Normal
Df Sum of Sq RSS AIC
- Electrical_SBrkr 1 0.00020 28.267 -5525.8
- ExterCond 1 0.00142 28.268 -5525.7
- YrSold_X2010 1 0.00507 28.272 -5525.6
- New 1 0.00640 28.273 -5525.5
- Fence_MnPrv 1 0.01335 28.280 -5525.1
- MoSold_Feb_Aug_Sep_Nov_Dec 1 0.01339 28.280 -5525.1
- BsmtCond 1 0.02193 28.289 -5524.7
- LandSlope_Mod_Sev 1 0.03222 28.299 -5524.2
- MSSubClass_X.150k 1 0.03929 28.306 -5523.8
- `2ndFlrSF` 1 0.05128 28.318 -5523.2
- Alley_None 1 0.06154 28.328 -5522.6
- LandContour_Lvl 1 0.06199 28.329 -5522.6
- RoofStyle_Gam_Gable 1 0.08334 28.350 -5521.5
- `1stFlrSF` 1 0.10349 28.370 -5520.5
<none> 28.267 -5518.5
- PavedDrive_Y 1 0.17313 28.440 -5516.9
- Neighborhood_X150k.200k 1 0.22541 28.492 -5514.2
- BldgType_fmCon_Dup_Twnhs 1 0.22971 28.497 -5514.0
- LotArea 1 0.25492 28.522 -5512.7
- GarageQual 1 0.27691 28.544 -5511.6
- CentralAir_Y 1 0.28673 28.553 -5511.1
- SaleCondition_Normal 1 0.33788 28.605 -5508.5
- RemodAdd 1 0.33860 28.605 -5508.4
- HeatingQC 1 0.40985 28.677 -5504.8
- MSZoning_RL 1 0.52665 28.794 -5498.9
- FireplaceQu 1 0.63457 28.901 -5493.4
- Age 1 0.97330 29.240 -5476.4
- Functional_Typ 1 1.04383 29.311 -5472.9
- GarageCars 1 1.58183 29.849 -5446.3
- GrLivArea 1 1.93397 30.201 -5429.2
- OverallCond 1 1.98051 30.247 -5426.9
- BsmtQual 1 2.99256 31.259 -5378.9
- Neighborhood_X.200k 1 3.15151 31.418 -5371.5
Step: AIC=-5525.81
SalePrice ~ LotArea + OverallCond + ExterCond + BsmtQual + BsmtCond +
HeatingQC + `1stFlrSF` + `2ndFlrSF` + GrLivArea + FireplaceQu +
GarageCars + GarageQual + RemodAdd + New + Age + MSSubClass_X.150k +
MSZoning_RL + Alley_None + LandContour_Lvl + LandSlope_Mod_Sev +
Neighborhood_X150k.200k + Neighborhood_X.200k + BldgType_fmCon_Dup_Twnhs +
RoofStyle_Gam_Gable + CentralAir_Y + Functional_Typ + PavedDrive_Y +
Fence_MnPrv + MoSold_Feb_Aug_Sep_Nov_Dec + YrSold_X2010 +
SaleCondition_Normal
Df Sum of Sq RSS AIC
- ExterCond 1 0.00134 28.268 -5533.0
- YrSold_X2010 1 0.00516 28.272 -5532.8
- New 1 0.00633 28.273 -5532.8
- Fence_MnPrv 1 0.01320 28.280 -5532.4
- MoSold_Feb_Aug_Sep_Nov_Dec 1 0.01345 28.280 -5532.4
- BsmtCond 1 0.02218 28.289 -5531.9
- LandSlope_Mod_Sev 1 0.03235 28.299 -5531.4
- MSSubClass_X.150k 1 0.03930 28.306 -5531.1
- `2ndFlrSF` 1 0.05148 28.319 -5530.4
- Alley_None 1 0.06159 28.329 -5529.9
- LandContour_Lvl 1 0.06220 28.329 -5529.9
- RoofStyle_Gam_Gable 1 0.08345 28.351 -5528.8
- `1stFlrSF` 1 0.10329 28.370 -5527.8
<none> 28.267 -5525.8
- PavedDrive_Y 1 0.17396 28.441 -5524.1
- Neighborhood_X150k.200k 1 0.22748 28.494 -5521.4
- BldgType_fmCon_Dup_Twnhs 1 0.22963 28.497 -5521.3
- LotArea 1 0.25472 28.522 -5520.0
- GarageQual 1 0.27719 28.544 -5518.8
- CentralAir_Y 1 0.29888 28.566 -5517.7
- SaleCondition_Normal 1 0.33835 28.605 -5515.7
- RemodAdd 1 0.34266 28.610 -5515.5
- HeatingQC 1 0.40982 28.677 -5512.1
- MSZoning_RL 1 0.52933 28.796 -5506.0
- FireplaceQu 1 0.63685 28.904 -5500.6
- Age 1 0.99837 29.265 -5482.4
- Functional_Typ 1 1.04408 29.311 -5480.1
- GarageCars 1 1.58615 29.853 -5453.4
- GrLivArea 1 1.93726 30.204 -5436.3
- OverallCond 1 1.99445 30.262 -5433.6
- BsmtQual 1 3.00017 31.267 -5385.8
- Neighborhood_X.200k 1 3.15292 31.420 -5378.7
Step: AIC=-5533.02
SalePrice ~ LotArea + OverallCond + BsmtQual + BsmtCond + HeatingQC +
`1stFlrSF` + `2ndFlrSF` + GrLivArea + FireplaceQu + GarageCars +
GarageQual + RemodAdd + New + Age + MSSubClass_X.150k + MSZoning_RL +
Alley_None + LandContour_Lvl + LandSlope_Mod_Sev + Neighborhood_X150k.200k +
Neighborhood_X.200k + BldgType_fmCon_Dup_Twnhs + RoofStyle_Gam_Gable +
CentralAir_Y + Functional_Typ + PavedDrive_Y + Fence_MnPrv +
MoSold_Feb_Aug_Sep_Nov_Dec + YrSold_X2010 + SaleCondition_Normal
Df Sum of Sq RSS AIC
- YrSold_X2010 1 0.00516 28.273 -5540.0
- New 1 0.00663 28.275 -5540.0
- Fence_MnPrv 1 0.01288 28.281 -5539.6
- MoSold_Feb_Aug_Sep_Nov_Dec 1 0.01359 28.282 -5539.6
- BsmtCond 1 0.02115 28.290 -5539.2
- LandSlope_Mod_Sev 1 0.03294 28.301 -5538.6
- MSSubClass_X.150k 1 0.03906 28.307 -5538.3
- `2ndFlrSF` 1 0.05255 28.321 -5537.6
- Alley_None 1 0.06114 28.329 -5537.2
- LandContour_Lvl 1 0.06223 28.331 -5537.1
- RoofStyle_Gam_Gable 1 0.08356 28.352 -5536.0
- `1stFlrSF` 1 0.10220 28.371 -5535.0
<none> 28.268 -5533.0
- PavedDrive_Y 1 0.17270 28.441 -5531.4
- Neighborhood_X150k.200k 1 0.22782 28.496 -5528.6
- BldgType_fmCon_Dup_Twnhs 1 0.22868 28.497 -5528.5
- LotArea 1 0.25391 28.522 -5527.3
- GarageQual 1 0.27594 28.544 -5526.1
- CentralAir_Y 1 0.29946 28.568 -5524.9
- SaleCondition_Normal 1 0.33894 28.607 -5522.9
- RemodAdd 1 0.34162 28.610 -5522.8
- HeatingQC 1 0.40887 28.677 -5519.3
- MSZoning_RL 1 0.53406 28.802 -5513.0
- FireplaceQu 1 0.63748 28.906 -5507.8
- Age 1 0.99715 29.265 -5489.7
- Functional_Typ 1 1.04323 29.312 -5487.4
- GarageCars 1 1.58557 29.854 -5460.6
- GrLivArea 1 1.94693 30.215 -5443.1
- OverallCond 1 2.17067 30.439 -5432.3
- BsmtQual 1 3.01602 31.284 -5392.3
- Neighborhood_X.200k 1 3.15315 31.422 -5385.9
Step: AIC=-5540.04
SalePrice ~ LotArea + OverallCond + BsmtQual + BsmtCond + HeatingQC +
`1stFlrSF` + `2ndFlrSF` + GrLivArea + FireplaceQu + GarageCars +
GarageQual + RemodAdd + New + Age + MSSubClass_X.150k + MSZoning_RL +
Alley_None + LandContour_Lvl + LandSlope_Mod_Sev + Neighborhood_X150k.200k +
Neighborhood_X.200k + BldgType_fmCon_Dup_Twnhs + RoofStyle_Gam_Gable +
CentralAir_Y + Functional_Typ + PavedDrive_Y + Fence_MnPrv +
MoSold_Feb_Aug_Sep_Nov_Dec + SaleCondition_Normal
Df Sum of Sq RSS AIC
- New 1 0.00637 28.280 -5547.0
- Fence_MnPrv 1 0.01263 28.286 -5546.7
- MoSold_Feb_Aug_Sep_Nov_Dec 1 0.01625 28.290 -5546.5
- BsmtCond 1 0.02068 28.294 -5546.3
- LandSlope_Mod_Sev 1 0.03332 28.307 -5545.6
- MSSubClass_X.150k 1 0.03920 28.313 -5545.3
- `2ndFlrSF` 1 0.05331 28.327 -5544.6
- Alley_None 1 0.06237 28.336 -5544.1
- LandContour_Lvl 1 0.06342 28.337 -5544.1
- RoofStyle_Gam_Gable 1 0.08270 28.356 -5543.1
- `1stFlrSF` 1 0.10160 28.375 -5542.1
<none> 28.273 -5540.0
- PavedDrive_Y 1 0.17423 28.448 -5538.4
- Neighborhood_X150k.200k 1 0.22636 28.500 -5535.7
- BldgType_fmCon_Dup_Twnhs 1 0.22813 28.502 -5535.6
- LotArea 1 0.25428 28.528 -5534.3
- GarageQual 1 0.27578 28.549 -5533.2
- CentralAir_Y 1 0.30002 28.573 -5531.9
- SaleCondition_Normal 1 0.34136 28.615 -5529.8
- RemodAdd 1 0.34178 28.615 -5529.8
- HeatingQC 1 0.40673 28.680 -5526.5
- MSZoning_RL 1 0.53554 28.809 -5519.9
- FireplaceQu 1 0.63340 28.907 -5515.0
- Age 1 0.99472 29.268 -5496.8
- Functional_Typ 1 1.04195 29.316 -5494.5
- GarageCars 1 1.58222 29.856 -5467.8
- GrLivArea 1 1.95301 30.227 -5449.8
- OverallCond 1 2.17522 30.449 -5439.1
- BsmtQual 1 3.02848 31.302 -5398.8
- Neighborhood_X.200k 1 3.15576 31.429 -5392.8
Step: AIC=-5547
SalePrice ~ LotArea + OverallCond + BsmtQual + BsmtCond + HeatingQC +
`1stFlrSF` + `2ndFlrSF` + GrLivArea + FireplaceQu + GarageCars +
GarageQual + RemodAdd + Age + MSSubClass_X.150k + MSZoning_RL +
Alley_None + LandContour_Lvl + LandSlope_Mod_Sev + Neighborhood_X150k.200k +
Neighborhood_X.200k + BldgType_fmCon_Dup_Twnhs + RoofStyle_Gam_Gable +
CentralAir_Y + Functional_Typ + PavedDrive_Y + Fence_MnPrv +
MoSold_Feb_Aug_Sep_Nov_Dec + SaleCondition_Normal
Df Sum of Sq RSS AIC
- Fence_MnPrv 1 0.01289 28.293 -5553.6
- MoSold_Feb_Aug_Sep_Nov_Dec 1 0.01468 28.295 -5553.5
- BsmtCond 1 0.02236 28.302 -5553.1
- LandSlope_Mod_Sev 1 0.03181 28.312 -5552.6
- MSSubClass_X.150k 1 0.03872 28.319 -5552.3
- `2ndFlrSF` 1 0.05226 28.332 -5551.6
- LandContour_Lvl 1 0.06144 28.341 -5551.1
- Alley_None 1 0.06193 28.342 -5551.1
- RoofStyle_Gam_Gable 1 0.08418 28.364 -5549.9
- `1stFlrSF` 1 0.10241 28.382 -5549.0
<none> 28.280 -5547.0
- PavedDrive_Y 1 0.17320 28.453 -5545.4
- Neighborhood_X150k.200k 1 0.22173 28.502 -5542.9
- BldgType_fmCon_Dup_Twnhs 1 0.23447 28.514 -5542.2
- LotArea 1 0.26073 28.541 -5540.9
- GarageQual 1 0.27658 28.556 -5540.1
- CentralAir_Y 1 0.29528 28.575 -5539.1
- SaleCondition_Normal 1 0.35858 28.639 -5535.9
- RemodAdd 1 0.40182 28.682 -5533.7
- HeatingQC 1 0.40539 28.685 -5533.5
- MSZoning_RL 1 0.53381 28.814 -5527.0
- FireplaceQu 1 0.63335 28.913 -5521.9
- Functional_Typ 1 1.03985 29.320 -5501.6
- Age 1 1.12845 29.408 -5497.2
- GarageCars 1 1.57672 29.857 -5475.1
- GrLivArea 1 1.94981 30.230 -5456.9
- OverallCond 1 2.17327 30.453 -5446.2
- BsmtQual 1 3.02553 31.305 -5405.9
- Neighborhood_X.200k 1 3.15298 31.433 -5400.0
Step: AIC=-5553.62
SalePrice ~ LotArea + OverallCond + BsmtQual + BsmtCond + HeatingQC +
`1stFlrSF` + `2ndFlrSF` + GrLivArea + FireplaceQu + GarageCars +
GarageQual + RemodAdd + Age + MSSubClass_X.150k + MSZoning_RL +
Alley_None + LandContour_Lvl + LandSlope_Mod_Sev + Neighborhood_X150k.200k +
Neighborhood_X.200k + BldgType_fmCon_Dup_Twnhs + RoofStyle_Gam_Gable +
CentralAir_Y + Functional_Typ + PavedDrive_Y + MoSold_Feb_Aug_Sep_Nov_Dec +
SaleCondition_Normal
Df Sum of Sq RSS AIC
- MoSold_Feb_Aug_Sep_Nov_Dec 1 0.0155 28.308 -5560.1
- BsmtCond 1 0.0211 28.314 -5559.8
- LandSlope_Mod_Sev 1 0.0321 28.325 -5559.2
- MSSubClass_X.150k 1 0.0391 28.332 -5558.9
- `2ndFlrSF` 1 0.0529 28.346 -5558.2
- LandContour_Lvl 1 0.0588 28.352 -5557.9
- Alley_None 1 0.0618 28.355 -5557.7
- RoofStyle_Gam_Gable 1 0.0829 28.376 -5556.6
- `1stFlrSF` 1 0.1007 28.393 -5555.7
<none> 28.293 -5553.6
- PavedDrive_Y 1 0.1715 28.464 -5552.1
- Neighborhood_X150k.200k 1 0.2220 28.515 -5549.5
- BldgType_fmCon_Dup_Twnhs 1 0.2257 28.518 -5549.3
- LotArea 1 0.2610 28.554 -5547.5
- GarageQual 1 0.2770 28.570 -5546.7
- CentralAir_Y 1 0.2904 28.583 -5546.0
- SaleCondition_Normal 1 0.3647 28.657 -5542.2
- RemodAdd 1 0.4094 28.702 -5539.9
- HeatingQC 1 0.4264 28.719 -5539.1
- MSZoning_RL 1 0.5292 28.822 -5533.8
- FireplaceQu 1 0.6418 28.934 -5528.2
- Functional_Typ 1 1.0614 29.354 -5507.1
- Age 1 1.1544 29.447 -5502.5
- GarageCars 1 1.5925 29.885 -5481.0
- GrLivArea 1 1.9596 30.252 -5463.1
- OverallCond 1 2.1604 30.453 -5453.5
- BsmtQual 1 3.0328 31.326 -5412.2
- Neighborhood_X.200k 1 3.1684 31.461 -5405.9
Step: AIC=-5560.11
SalePrice ~ LotArea + OverallCond + BsmtQual + BsmtCond + HeatingQC +
`1stFlrSF` + `2ndFlrSF` + GrLivArea + FireplaceQu + GarageCars +
GarageQual + RemodAdd + Age + MSSubClass_X.150k + MSZoning_RL +
Alley_None + LandContour_Lvl + LandSlope_Mod_Sev + Neighborhood_X150k.200k +
Neighborhood_X.200k + BldgType_fmCon_Dup_Twnhs + RoofStyle_Gam_Gable +
CentralAir_Y + Functional_Typ + PavedDrive_Y + SaleCondition_Normal
Df Sum of Sq RSS AIC
- BsmtCond 1 0.01902 28.327 -5566.4
- LandSlope_Mod_Sev 1 0.03193 28.340 -5565.7
- MSSubClass_X.150k 1 0.03729 28.346 -5565.5
- `2ndFlrSF` 1 0.05455 28.363 -5564.6
- LandContour_Lvl 1 0.06185 28.370 -5564.2
- Alley_None 1 0.06240 28.371 -5564.2
- RoofStyle_Gam_Gable 1 0.08440 28.393 -5563.0
- `1stFlrSF` 1 0.09833 28.407 -5562.3
<none> 28.308 -5560.1
- PavedDrive_Y 1 0.17043 28.479 -5558.6
- Neighborhood_X150k.200k 1 0.21608 28.524 -5556.3
- BldgType_fmCon_Dup_Twnhs 1 0.22927 28.538 -5555.6
- LotArea 1 0.26051 28.569 -5554.0
- GarageQual 1 0.27898 28.587 -5553.1
- CentralAir_Y 1 0.29039 28.599 -5552.5
- SaleCondition_Normal 1 0.37453 28.683 -5548.2
- RemodAdd 1 0.40543 28.714 -5546.6
- HeatingQC 1 0.43295 28.741 -5545.2
- MSZoning_RL 1 0.53667 28.845 -5540.0
- FireplaceQu 1 0.63838 28.947 -5534.8
- Functional_Typ 1 1.06669 29.375 -5513.4
- Age 1 1.14907 29.457 -5509.3
- GarageCars 1 1.60586 29.914 -5486.8
- GrLivArea 1 1.96752 30.276 -5469.3
- OverallCond 1 2.17026 30.479 -5459.5
- BsmtQual 1 3.03044 31.339 -5418.9
- Neighborhood_X.200k 1 3.15411 31.462 -5413.2
Step: AIC=-5566.41
SalePrice ~ LotArea + OverallCond + BsmtQual + HeatingQC + `1stFlrSF` +
`2ndFlrSF` + GrLivArea + FireplaceQu + GarageCars + GarageQual +
RemodAdd + Age + MSSubClass_X.150k + MSZoning_RL + Alley_None +
LandContour_Lvl + LandSlope_Mod_Sev + Neighborhood_X150k.200k +
Neighborhood_X.200k + BldgType_fmCon_Dup_Twnhs + RoofStyle_Gam_Gable +
CentralAir_Y + Functional_Typ + PavedDrive_Y + SaleCondition_Normal
Df Sum of Sq RSS AIC
- MSSubClass_X.150k 1 0.0358 28.363 -5571.9
- LandSlope_Mod_Sev 1 0.0363 28.364 -5571.8
- `2ndFlrSF` 1 0.0554 28.383 -5570.8
- Alley_None 1 0.0598 28.387 -5570.6
- LandContour_Lvl 1 0.0661 28.393 -5570.3
- RoofStyle_Gam_Gable 1 0.0829 28.410 -5569.4
- `1stFlrSF` 1 0.0963 28.424 -5568.7
<none> 28.327 -5566.4
- PavedDrive_Y 1 0.1769 28.504 -5564.6
- Neighborhood_X150k.200k 1 0.2105 28.538 -5562.9
- BldgType_fmCon_Dup_Twnhs 1 0.2344 28.562 -5561.7
- LotArea 1 0.2648 28.592 -5560.1
- GarageQual 1 0.2789 28.606 -5559.4
- CentralAir_Y 1 0.3085 28.636 -5557.9
- SaleCondition_Normal 1 0.3654 28.693 -5555.0
- RemodAdd 1 0.4118 28.739 -5552.6
- HeatingQC 1 0.4320 28.759 -5551.6
- MSZoning_RL 1 0.5416 28.869 -5546.1
- FireplaceQu 1 0.6371 28.964 -5541.2
- Functional_Typ 1 1.0810 29.408 -5519.0
- Age 1 1.1342 29.461 -5516.4
- GarageCars 1 1.6032 29.930 -5493.3
- GrLivArea 1 1.9723 30.300 -5475.4
- OverallCond 1 2.3591 30.686 -5456.9
- Neighborhood_X.200k 1 3.1352 31.462 -5420.4
- BsmtQual 1 3.6434 31.971 -5397.0
Step: AIC=-5571.86
SalePrice ~ LotArea + OverallCond + BsmtQual + HeatingQC + `1stFlrSF` +
`2ndFlrSF` + GrLivArea + FireplaceQu + GarageCars + GarageQual +
RemodAdd + Age + MSZoning_RL + Alley_None + LandContour_Lvl +
LandSlope_Mod_Sev + Neighborhood_X150k.200k + Neighborhood_X.200k +
BldgType_fmCon_Dup_Twnhs + RoofStyle_Gam_Gable + CentralAir_Y +
Functional_Typ + PavedDrive_Y + SaleCondition_Normal
Df Sum of Sq RSS AIC
- LandSlope_Mod_Sev 1 0.0348 28.398 -5577.4
- Alley_None 1 0.0675 28.431 -5575.7
- LandContour_Lvl 1 0.0690 28.432 -5575.6
- `1stFlrSF` 1 0.0761 28.439 -5575.2
- `2ndFlrSF` 1 0.0830 28.446 -5574.9
- RoofStyle_Gam_Gable 1 0.0917 28.455 -5574.4
<none> 28.363 -5571.9
- PavedDrive_Y 1 0.1834 28.547 -5569.7
- Neighborhood_X150k.200k 1 0.2250 28.588 -5567.6
- BldgType_fmCon_Dup_Twnhs 1 0.2544 28.617 -5566.1
- GarageQual 1 0.2768 28.640 -5565.0
- CentralAir_Y 1 0.3020 28.665 -5563.7
- LotArea 1 0.3418 28.705 -5561.7
- SaleCondition_Normal 1 0.3601 28.723 -5560.7
- HeatingQC 1 0.4229 28.786 -5557.5
- RemodAdd 1 0.4519 28.815 -5556.1
- MSZoning_RL 1 0.5966 28.960 -5548.8
- FireplaceQu 1 0.6122 28.975 -5548.0
- Functional_Typ 1 1.1342 29.497 -5521.9
- Age 1 1.1918 29.555 -5519.0
- GarageCars 1 1.6328 29.996 -5497.4
- GrLivArea 1 2.2809 30.644 -5466.2
- OverallCond 1 2.3719 30.735 -5461.9
- Neighborhood_X.200k 1 3.1067 31.470 -5427.4
- BsmtQual 1 3.6077 31.971 -5404.3
Step: AIC=-5577.35
SalePrice ~ LotArea + OverallCond + BsmtQual + HeatingQC + `1stFlrSF` +
`2ndFlrSF` + GrLivArea + FireplaceQu + GarageCars + GarageQual +
RemodAdd + Age + MSZoning_RL + Alley_None + LandContour_Lvl +
Neighborhood_X150k.200k + Neighborhood_X.200k + BldgType_fmCon_Dup_Twnhs +
RoofStyle_Gam_Gable + CentralAir_Y + Functional_Typ + PavedDrive_Y +
SaleCondition_Normal
Df Sum of Sq RSS AIC
- LandContour_Lvl 1 0.0383 28.436 -5582.7
- Alley_None 1 0.0709 28.469 -5581.0
- `1stFlrSF` 1 0.0757 28.474 -5580.8
- `2ndFlrSF` 1 0.0828 28.481 -5580.4
- RoofStyle_Gam_Gable 1 0.0917 28.490 -5579.9
<none> 28.398 -5577.4
- PavedDrive_Y 1 0.1894 28.587 -5574.9
- Neighborhood_X150k.200k 1 0.2311 28.629 -5572.8
- BldgType_fmCon_Dup_Twnhs 1 0.2540 28.652 -5571.6
- GarageQual 1 0.2787 28.677 -5570.4
- CentralAir_Y 1 0.3083 28.706 -5568.9
- SaleCondition_Normal 1 0.3638 28.762 -5566.1
- LotArea 1 0.3738 28.772 -5565.5
- HeatingQC 1 0.4211 28.819 -5563.1
- RemodAdd 1 0.4524 28.850 -5561.6
- MSZoning_RL 1 0.5869 28.985 -5554.8
- FireplaceQu 1 0.6097 29.008 -5553.6
- Functional_Typ 1 1.1050 29.503 -5528.9
- Age 1 1.1746 29.572 -5525.5
- GarageCars 1 1.6234 30.021 -5503.5
- GrLivArea 1 2.2737 30.672 -5472.2
- OverallCond 1 2.3779 30.776 -5467.2
- Neighborhood_X.200k 1 3.1839 31.582 -5429.5
- BsmtQual 1 3.6144 32.012 -5409.7
Step: AIC=-5582.67
SalePrice ~ LotArea + OverallCond + BsmtQual + HeatingQC + `1stFlrSF` +
`2ndFlrSF` + GrLivArea + FireplaceQu + GarageCars + GarageQual +
RemodAdd + Age + MSZoning_RL + Alley_None + Neighborhood_X150k.200k +
Neighborhood_X.200k + BldgType_fmCon_Dup_Twnhs + RoofStyle_Gam_Gable +
CentralAir_Y + Functional_Typ + PavedDrive_Y + SaleCondition_Normal
Df Sum of Sq RSS AIC
- `1stFlrSF` 1 0.0726 28.509 -5586.2
- Alley_None 1 0.0751 28.511 -5586.1
- `2ndFlrSF` 1 0.0868 28.523 -5585.5
- RoofStyle_Gam_Gable 1 0.0906 28.527 -5585.3
<none> 28.436 -5582.7
- PavedDrive_Y 1 0.2008 28.637 -5579.7
- Neighborhood_X150k.200k 1 0.2206 28.657 -5578.7
- BldgType_fmCon_Dup_Twnhs 1 0.2683 28.705 -5576.2
- GarageQual 1 0.2833 28.720 -5575.5
- CentralAir_Y 1 0.3107 28.747 -5574.1
- LotArea 1 0.3426 28.779 -5572.5
- SaleCondition_Normal 1 0.3732 28.809 -5570.9
- HeatingQC 1 0.4342 28.870 -5567.8
- RemodAdd 1 0.4627 28.899 -5566.4
- FireplaceQu 1 0.5972 29.033 -5559.6
- MSZoning_RL 1 0.6179 29.054 -5558.6
- Functional_Typ 1 1.1186 29.555 -5533.6
- Age 1 1.1928 29.629 -5530.0
- GarageCars 1 1.6344 30.071 -5508.4
- GrLivArea 1 2.3120 30.748 -5475.8
- OverallCond 1 2.3643 30.801 -5473.3
- Neighborhood_X.200k 1 3.1516 31.588 -5436.5
- BsmtQual 1 3.6052 32.041 -5415.7
Step: AIC=-5586.23
SalePrice ~ LotArea + OverallCond + BsmtQual + HeatingQC + `2ndFlrSF` +
GrLivArea + FireplaceQu + GarageCars + GarageQual + RemodAdd +
Age + MSZoning_RL + Alley_None + Neighborhood_X150k.200k +
Neighborhood_X.200k + BldgType_fmCon_Dup_Twnhs + RoofStyle_Gam_Gable +
CentralAir_Y + Functional_Typ + PavedDrive_Y + SaleCondition_Normal
Df Sum of Sq RSS AIC
- Alley_None 1 0.0764 28.585 -5589.6
- RoofStyle_Gam_Gable 1 0.0887 28.598 -5589.0
<none> 28.509 -5586.2
- PavedDrive_Y 1 0.1960 28.705 -5583.5
- Neighborhood_X150k.200k 1 0.1999 28.709 -5583.3
- BldgType_fmCon_Dup_Twnhs 1 0.2748 28.784 -5579.5
- GarageQual 1 0.2902 28.799 -5578.7
- CentralAir_Y 1 0.3214 28.830 -5577.1
- SaleCondition_Normal 1 0.3624 28.871 -5575.1
- LotArea 1 0.3915 28.900 -5573.6
- HeatingQC 1 0.4295 28.938 -5571.7
- RemodAdd 1 0.4403 28.949 -5571.1
- FireplaceQu 1 0.6183 29.127 -5562.2
- MSZoning_RL 1 0.6194 29.128 -5562.1
- Functional_Typ 1 1.0849 29.594 -5539.0
- Age 1 1.1662 29.675 -5535.0
- GarageCars 1 1.6322 30.141 -5512.2
- `2ndFlrSF` 1 2.0406 30.549 -5492.6
- OverallCond 1 2.4071 30.916 -5475.2
- Neighborhood_X.200k 1 3.1062 31.615 -5442.5
- BsmtQual 1 3.6271 32.136 -5418.7
- GrLivArea 1 17.3706 45.879 -4898.8
Step: AIC=-5589.61
SalePrice ~ LotArea + OverallCond + BsmtQual + HeatingQC + `2ndFlrSF` +
GrLivArea + FireplaceQu + GarageCars + GarageQual + RemodAdd +
Age + MSZoning_RL + Neighborhood_X150k.200k + Neighborhood_X.200k +
BldgType_fmCon_Dup_Twnhs + RoofStyle_Gam_Gable + CentralAir_Y +
Functional_Typ + PavedDrive_Y + SaleCondition_Normal
Df Sum of Sq RSS AIC
- RoofStyle_Gam_Gable 1 0.0959 28.681 -5592.0
<none> 28.585 -5589.6
- Neighborhood_X150k.200k 1 0.2151 28.800 -5586.0
- PavedDrive_Y 1 0.2319 28.817 -5585.1
- BldgType_fmCon_Dup_Twnhs 1 0.2686 28.854 -5583.2
- GarageQual 1 0.2872 28.872 -5582.3
- CentralAir_Y 1 0.3544 28.940 -5578.9
- SaleCondition_Normal 1 0.3742 28.959 -5577.9
- HeatingQC 1 0.4151 29.000 -5575.8
- LotArea 1 0.4220 29.007 -5575.5
- RemodAdd 1 0.4618 29.047 -5573.5
- FireplaceQu 1 0.6560 29.241 -5563.8
- MSZoning_RL 1 0.6669 29.252 -5563.2
- Functional_Typ 1 1.0892 29.674 -5542.3
- Age 1 1.1631 29.748 -5538.7
- GarageCars 1 1.6045 30.190 -5517.2
- `2ndFlrSF` 1 2.1151 30.700 -5492.7
- OverallCond 1 2.3804 30.966 -5480.1
- Neighborhood_X.200k 1 3.0737 31.659 -5447.8
- BsmtQual 1 3.6742 32.259 -5420.4
- GrLivArea 1 17.3370 45.922 -4904.8
Step: AIC=-5592.01
SalePrice ~ LotArea + OverallCond + BsmtQual + HeatingQC + `2ndFlrSF` +
GrLivArea + FireplaceQu + GarageCars + GarageQual + RemodAdd +
Age + MSZoning_RL + Neighborhood_X150k.200k + Neighborhood_X.200k +
BldgType_fmCon_Dup_Twnhs + CentralAir_Y + Functional_Typ +
PavedDrive_Y + SaleCondition_Normal
Df Sum of Sq RSS AIC
<none> 28.681 -5592.0
- Neighborhood_X150k.200k 1 0.1852 28.866 -5589.9
- PavedDrive_Y 1 0.2269 28.908 -5587.8
- BldgType_fmCon_Dup_Twnhs 1 0.2758 28.957 -5585.3
- GarageQual 1 0.2927 28.974 -5584.5
- CentralAir_Y 1 0.3404 29.021 -5582.1
- SaleCondition_Normal 1 0.3796 29.061 -5580.1
- HeatingQC 1 0.4104 29.092 -5578.6
- LotArea 1 0.4502 29.131 -5576.6
- RemodAdd 1 0.4972 29.178 -5574.2
- MSZoning_RL 1 0.6623 29.343 -5566.0
- FireplaceQu 1 0.6813 29.362 -5565.0
- Functional_Typ 1 1.0742 29.755 -5545.6
- Age 1 1.1514 29.832 -5541.8
- GarageCars 1 1.6213 30.302 -5519.0
- OverallCond 1 2.4337 31.115 -5480.4
- `2ndFlrSF` 1 2.4685 31.150 -5478.8
- Neighborhood_X.200k 1 3.0523 31.733 -5451.6
- BsmtQual 1 3.7266 32.408 -5420.9
- GrLivArea 1 18.2085 46.890 -4881.6
Resumimos el modelo tras su procesado con summary():
summary(modBIC)
Call:
lm(formula = SalePrice ~ LotArea + OverallCond + BsmtQual + HeatingQC +
`2ndFlrSF` + GrLivArea + FireplaceQu + GarageCars + GarageQual +
RemodAdd + Age + MSZoning_RL + Neighborhood_X150k.200k +
Neighborhood_X.200k + BldgType_fmCon_Dup_Twnhs + CentralAir_Y +
Functional_Typ + PavedDrive_Y + SaleCondition_Normal, data = house_prep)
Residuals:
Min 1Q Median 3Q Max
-1.13018 -0.07507 0.00400 0.08492 0.50303
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 11.663334 0.026693 436.937 < 2e-16 ***
LotArea 0.022914 0.004820 4.754 2.19e-06 ***
OverallCond 0.047640 0.004310 11.054 < 2e-16 ***
BsmtQual 0.074346 0.005435 13.679 < 2e-16 ***
HeatingQC 0.021426 0.004720 4.539 6.11e-06 ***
`2ndFlrSF` -0.057006 0.005121 -11.133 < 2e-16 ***
GrLivArea 0.201236 0.006656 30.236 < 2e-16 ***
FireplaceQu 0.026632 0.004554 5.849 6.12e-09 ***
GarageCars 0.049482 0.005485 9.022 < 2e-16 ***
GarageQual 0.016549 0.004317 3.834 0.000132 ***
RemodAdd -0.021793 0.004362 -4.996 6.56e-07 ***
Age -0.042503 0.005590 -7.603 5.18e-14 ***
MSZoning_RL 0.070250 0.012182 5.766 9.89e-09 ***
Neighborhood_X150k.200k 0.033704 0.011054 3.049 0.002337 **
Neighborhood_X.200k 0.157659 0.012736 12.379 < 2e-16 ***
BldgType_fmCon_Dup_Twnhs -0.054103 0.014539 -3.721 0.000206 ***
CentralAir_Y 0.072612 0.017565 4.134 3.77e-05 ***
Functional_Typ 0.114182 0.015548 7.344 3.46e-13 ***
PavedDrive_Y 0.052063 0.015426 3.375 0.000758 ***
SaleCondition_Normal 0.044835 0.010270 4.365 1.36e-05 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.1411 on 1440 degrees of freedom
Multiple R-squared: 0.8768, Adjusted R-squared: 0.8752
F-statistic: 539.4 on 19 and 1440 DF, p-value: < 2.2e-16
Como se puede observar, ha eliminado una gran cantidad de variables. Pasamos ahora a reescribir la receta únicamente con las variables seleccionadas:
house_bic <-
house_train |> mutate(MSZoning = MSZoning == "RL")
house_bic <-
house_bic |> mutate(Neighborhood = Neighborhood == "X.200k" )
house_bic <-
house_bic |> mutate(BldgType = BldgType == "fmCon_Dup_Twnhs" )
house_bic <-
house_bic |> mutate(CentralAir = CentralAir == "Y")
house_bic <-
house_bic |> mutate(Functional = Functional == "Typ")
house_bic <-
house_bic |> mutate(PavedDrive = PavedDrive == "Y")
house_bic <-
house_bic |> mutate(SaleCondition = SaleCondition == "Normal")
Volvemos a definir una nueva receta indicándole el conjunto donde tenemos validación y train, y enfrentando nuestra variable objetivo SalePrice a las variables seleccionadas por step_AIC(). Después, asignamos posibles roles, sujetos a modificación, que nos permitan diferenciar acciones entre las variables (sobre todo en la sección outliers).
# Receta
house_rec_multi_bic <-
# Fórmula y datos
recipe(data = house_bic, SalePrice ~ LotArea + OverallCond + BsmtQual +
HeatingQC + `2ndFlrSF` + GrLivArea + FireplaceQu + GarageCars +
GarageQual + RemodAdd + Age + MSZoning + Neighborhood + BldgType +
CentralAir + Functional + PavedDrive + SaleCondition)|>
# Roles
add_role(where(is.factor),
new_role = "cualitativa") |>
add_role(where(is.numeric),
new_role = "cuantitativa") |>
add_role(LotArea, `2ndFlrSF`,
new_role = "mediana") |>
add_role(all_nominal_predictors(),
new_role = "moda") |>
add_role(OverallCond,
new_role = "imputar media")
En este modelo no será necesario reagrupar. El criterio de información BIC ya ha seleccionado las categorías más apropiadas, si las hubiera, de cada variable.
En los métodos de regresión, la recategorización de las variables cuantitativas no es estrictamente necesaria. Por tanto, en este receta, del total de variables cuantitativas no vamos a modificar ninguna.
multiplot(box1, box2, box3, box4, box5, box6, box7, box8, box9, box10, cols = 2)
multiplot(box11, box12, box13, box14, box15, box16, box17, cols = 2)
Si observamos estos gráficos de cajas y bigotes, todas nuestras variables cuantitativas continuas son asimétricas, por lo que se detectarán los outliers y se imputarán los ausentes por la mediana. Para el caso de las semi-cuantitativas (KitchenQual, por ejemplo) se imputarán directamente los ausentes por uno de los dos métodos en función de la simetría o asimetría de la variable. Para el resto de variables cualitativas, imputamos directamente por la moda.
house_rec_multi_bic <-
house_rec_multi_bic |>
# Detección de outliers por la mediana y por la media
step_mutate(across(has_role("mediana"), function(x) { ifelse(abs(scores(x, type = "mad")) > 3 & !is.na(x), NA, x) })) |>
# Imputación de ausentes por la mediana, la media y la moda
step_impute_median(has_role("mediana")) |>
step_impute_mean(has_role("imputar_media")) |>
step_impute_mode(has_role("moda")) |>
step_impute_knn(all_numeric_predictors() , -has_role("mediana"), -has_role("imputar_media") )
A fin de que nuestro modelo cumpla los test de normalidad de los residuos, vamos a incluir en la receta la función step_YeoJohnson, que se emplea para eliminar la asimetría en los datos numéricos del modelo.
house_rec_multi_bic <-
house_rec_multi_bic |>
step_YeoJohnson(all_numeric_predictors())
Transformamos algunas variables a logarítmicas para mejorar su relación.
house_rec_multi_bic <-
house_rec_multi_bic |>
step_log(SalePrice, offset = 1)
Normalizamos nuestras variables por rango para que todas tengan el mismo peso, entre 0 y 1.
house_rec_multi_bic <-
house_rec_multi_bic |>
step_normalize(all_numeric_predictors())
Como esta receta es para un modelo de regresión multivariante, debemos dummyficar nuestras variables cualitativas. Para ello, tomamos todas las nominales, menos nuestra variable objetivo.
house_rec_multi_bic <-
house_rec_multi_bic |>
step_dummy(all_nominal(), -all_outcomes())
house_rec_multi_bic <-
house_rec_multi_bic |>
step_zv(all_predictors())
Aplicamos el filtro de correlación a nuestras variables numéricas.
house_rec_multi_bic <-
house_rec_multi_bic |>
step_corr(all_numeric_predictors(), threshold = 0.6)
Por último, horneamos nuestra receta para comprobar que todas nuestras nuevas variables recategorizadas se hayan creado correctamente.
# A tibble: 1,460 × 17
LotArea OverallCond BsmtQual HeatingQC `2ndFlrSF` GrLivArea
<dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
1 -0.141 -0.477 0.576 0.930 1.17 0.528
2 0.106 2.01 0.576 0.930 -0.871 -0.383
3 0.414 -0.477 0.576 0.930 1.17 0.659
4 0.0955 -0.477 -0.695 -0.330 1.15 0.541
5 0.877 -0.477 0.576 0.930 1.21 1.28
6 0.857 -0.477 0.576 0.930 1.09 -0.154
7 0.201 -0.477 2.13 0.930 -0.871 0.500
8 0.257 0.440 0.576 0.930 1.20 1.13
9 -0.761 -0.477 -0.695 -0.330 1.15 0.639
10 -0.391 0.440 -0.695 0.930 -0.871 -0.857
# … with 1,450 more rows, and 11 more variables: FireplaceQu <dbl>,
# GarageCars <dbl>, GarageQual <dbl>, RemodAdd <dbl>, Age <dbl>,
# MSZoning <lgl>, CentralAir <lgl>, Functional <lgl>,
# PavedDrive <lgl>, SaleCondition <lgl>, SalePrice <dbl>
# Creación del flujo
reg_wflow_bic <-
workflow() |>
add_model(reg_lineal) |>
add_recipe(house_rec_multi_bic)
set.seed(05492)
reg_fit_bic <-
reg_wflow_bic |> fit(data = house_train)
reg_fit_bic
══ Workflow [trained] ════════════════════════════════════════════════
Preprocessor: Recipe
Model: linear_reg()
── Preprocessor ──────────────────────────────────────────────────────
11 Recipe Steps
• step_mutate()
• step_impute_median()
• step_impute_mean()
• step_impute_mode()
• step_impute_knn()
• step_YeoJohnson()
• step_log()
• step_normalize()
• step_dummy()
• step_zv()
• ...
• and 1 more step.
── Model ─────────────────────────────────────────────────────────────
Call:
stats::lm(formula = ..y ~ ., data = data)
Coefficients:
(Intercept) LotArea OverallCond
11.8486103 0.0354377 0.0491536
BsmtQual HeatingQC `2ndFlrSF`
0.0673094 0.0183775 -0.0548430
GrLivArea FireplaceQu GarageCars
0.1960803 0.0267381 0.0432562
GarageQual RemodAdd Age
0.0157258 -0.0167712 -0.0378894
MSZoning_FV MSZoning_RH MSZoning_RM
0.0797640 0.0090258 0.0296974
Neighborhood_Blueste Neighborhood_BrDale Neighborhood_BrkSide
-0.1233748 -0.0445626 -0.1180328
Neighborhood_ClearCr Neighborhood_CollgCr Neighborhood_Crawfor
0.0032435 -0.0258185 0.0083973
Neighborhood_Edwards Neighborhood_Gilbert Neighborhood_IDOTRR
-0.1625234 -0.0772346 -0.2696271
Neighborhood_MeadowV Neighborhood_Mitchel Neighborhood_NAmes
-0.1879102 -0.0919997 -0.0930020
Neighborhood_NoRidge Neighborhood_NPkVill Neighborhood_NridgHt
0.1261509 -0.0229875 0.1328166
Neighborhood_NWAmes Neighborhood_OldTown Neighborhood_Sawyer
-0.1066500 -0.2251197 -0.1199513
Neighborhood_SawyerW Neighborhood_StoneBr Neighborhood_SWISU
-0.0555471 0.2130249 -0.1420178
Neighborhood_Timber Neighborhood_Veenker BldgType_X2fmCon
-0.0080007 0.0485150 -0.0215432
BldgType_Duplex BldgType_Twnhs BldgType_TwnhsE
-0.0860395 -0.0819172 -0.0496835
CentralAir_Y Functional_Maj2 Functional_Min1
0.0732900 -0.1285886 0.0054795
Functional_Min2 Functional_Mod Functional_Sev
0.0247964 0.0169903 -0.3563441
Functional_Typ PavedDrive_P PavedDrive_Y
0.1103448 0.0005071 0.0463707
SaleCondition_AdjLand SaleCondition_Alloca SaleCondition_Family
0.1887456 0.0061140 -0.0288674
SaleCondition_Normal
0.0385597
# Resumen del ajuste
tidy(reg_fit_bic)
# A tibble: 55 × 5
term estimate std.error statistic p.value
<chr> <dbl> <dbl> <dbl> <dbl>
1 (Intercept) 11.8 0.0533 222. 0
2 LotArea 0.0354 0.00609 5.82 7.32e- 9
3 OverallCond 0.0492 0.00455 10.8 3.37e- 26
4 BsmtQual 0.0673 0.00556 12.1 3.80e- 32
5 HeatingQC 0.0184 0.00487 3.77 1.67e- 4
6 `2ndFlrSF` -0.0548 0.00545 -10.1 4.38e- 23
7 GrLivArea 0.196 0.00684 28.7 1.11e-142
8 FireplaceQu 0.0267 0.00469 5.70 1.45e- 8
9 GarageCars 0.0433 0.00555 7.79 1.31e- 14
10 GarageQual 0.0157 0.00430 3.66 2.66e- 4
# … with 45 more rows
reg_fit_bic |> extract_fit_engine() |> performance()
# Indices of model performance
AIC | AICc | BIC | R2 | R2 (adj.) | RMSE | Sigma
---------------------------------------------------------------------
-1627.362 | -1622.812 | -1331.336 | 0.888 | 0.884 | 0.133 | 0.136
Nuestro \(R^2\) es igual a 0.888, por lo que el ratio de información explicada del modelo es del 89 %. Este modelo es superior al ratio de información explicada con el modelo multivariante anterior.
check_model(reg_fit_bic |> extract_fit_engine())
En el gráfico, la línea es bastante atendencial, se ve bastante monótona en todo su recorrido (excepto en sus extremos). Comprobémoslo analíticamente con un ANOVA entre residuos y predictora:
adjustment <-
reg_fit_bic |> extract_fit_engine()
lm(adjustment$residuals ~ adjustment$fitted.values) |> anova()
Analysis of Variance Table
Response: adjustment$residuals
Df Sum Sq Mean Sq F value Pr(>F)
adjustment$fitted.values 1 0.000 0.000000 0 1
Residuals 1458 25.971 0.017813
Analysis of Variance Table
Response: adjustment$residuals
Df Sum Sq Mean Sq F value Pr(>F)
I(adjustment$fitted.values^2) 1 0.0001 0.000061 0.0034 0.953411
adjustment$fitted.values 1 0.1195 0.119476 6.7338 0.009555
Residuals 1457 25.8512 0.017743
I(adjustment$fitted.values^2)
adjustment$fitted.values **
Residuals
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Como se puede observar, el p-valor de la prueba F de Fisher-Snedecor es igual a 1, por lo que, a un nivel de significación del 0.05, no podemos rechazar la hipótesis nula, esto es, la no existencia de una relación lineal entre las predictoras y sus residuos. Por tanto, ello implica que no parece existir tendencia lineal entre residuos y predictoras en nuestro modelo. Lo mismo ocurre en el segundo caso: tampoco parece existir tendencia cuadrática entre residuos y predictoras.
En el gráfico, la línea sigue una fase decreciente al principio y creciente al final, se ve poco monótona en todo su recorrido. Comprobémoslo analíticamente con un test de heterocedasticidad:
check_heteroscedasticity(adjustment)
Warning: Heteroscedasticity (non-constant error variance) detected (p < .001).
Según el test de heterocedasticidad, el supuesto de homocedasticidad no se cumple para nuestro modelo. Esto, muy a menudo, resulta poco realista, pues el supuesto implica que la variabilidad del error de la variable objetivo es la misma para cualquier nivel de nuestra predictora.
ggplot(
tibble("Observaciones" = 1:length(adjustment$residuals),
"Residuos" = adjustment$residuals),
aes(x = Observaciones, y = Residuos)) +
geom_point() +
geom_smooth(method = "lm", se = FALSE) +
theme_minimal()
En el gráfico, la recta de regresión se aprecia constante, en torno a 0, y los residuos se agrupan en torno a una banda también relativamente constante.
En el gráfico original de la diagnosis, prácticamente todas las observaciones están en la franja verde. Ninguna está en la franja roja.
En el gráfico, los percentiles empíricos de los residuos se acoplan bastante bien a la diagonal, excepto en los extremos. Comprobémoslo analíticamente con un contraste de normalidad:
library(olsrr)
ols_test_normality(adjustment$residuals)
-----------------------------------------------
Test Statistic pvalue
-----------------------------------------------
Shapiro-Wilk 0.939 0.0000
Kolmogorov-Smirnov 0.0616 0.0000
Cramer-von Mises 377.7789 0.0000
Anderson-Darling 10.0147 0.0000
-----------------------------------------------
Como se puede observar, los p-valor de los distintos test son iguales a 0, por lo que, a un nivel de significación del 0.05, podemos rechazar la hipótesis de normalidad. Por tanto, no podemos asumir a normalidad en nuestro modelo. Se han hecho pruebas con cambios en la receta añadiendo distintas transformaciones a nuestra objetivo (step_YeoJohnson(), step_BoxCox(), step_sqrt()), pero tampoco ha sido posible.
En el gráfico original de la diagnosis, todas las observaciones parecen estar contenidas entre las líneas punteadas al 0.5.
glance(reg_fit_bic)
# A tibble: 1 × 12
r.squared adj.r.squared sigma statistic p.value df logLik AIC
<dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
1 0.888 0.884 0.136 207. 0 54 870. -1627.
# … with 4 more variables: BIC <dbl>, deviance <dbl>,
# df.residual <int>, nobs <int>
Nuestro \(R^2\) es igual a 0.8884, por lo que el ratio de información explicada del modelo es del 89 %.
# Hacemos el split. Lo hacemos del 0.6 ya que hay pocos datos en train.
set.seed(05492)
split_house <-
initial_split(house_train, prop = 0.6)
# Predecimos en test
reg_fit_bic_test <-
reg_fit_bic |> last_fit(split = split_house)
# Evaluamos en test
reg_fit_bic_test |> collect_metrics()
# A tibble: 2 × 4
.metric .estimator .estimate .config
<chr> <chr> <dbl> <chr>
1 rmse standard 0.149 Preprocessor1_Model1
2 rsq standard 0.856 Preprocessor1_Model1
set.seed(05492)
# Errores en test
pred_test <-
reg_fit_bic_test |>
collect_predictions() |>
mutate(error = SalePrice - .pred)
pred_test
# A tibble: 584 × 6
id .pred .row SalePrice .config error
<chr> <dbl> <int> <dbl> <chr> <dbl>
1 train/test split 12.3 2 12.1 Preprocessor1_Model1 -0.152
2 train/test split 12.2 3 12.3 Preprocessor1_Model1 0.0711
3 train/test split 12.5 7 12.6 Preprocessor1_Model1 0.134
4 train/test split 12.3 8 12.2 Preprocessor1_Model1 -0.0565
5 train/test split 11.7 9 11.8 Preprocessor1_Model1 0.0988
6 train/test split 11.7 11 11.8 Preprocessor1_Model1 0.0742
7 train/test split 11.7 13 11.9 Preprocessor1_Model1 0.224
8 train/test split 11.9 15 12.0 Preprocessor1_Model1 0.106
9 train/test split 11.9 17 11.9 Preprocessor1_Model1 -0.0135
10 train/test split 11.7 20 11.8 Preprocessor1_Model1 0.146
# … with 574 more rows
En esta tabla se puede observar cómo, en nuestro modelo, las estimaciones se ajustan en gran medida a los valores reales de la objetivo SalePrice. Visualizaremos ahora estos valores a través de un gráfico.
g1 <- pred_test |>
ggplot(mapping = aes(x = .pred, y = SalePrice)) +
geom_point(color = "#56BCC2", alpha = 0.6, size = 4) +
geom_abline(intercept = 0, slope = 1, color = "#EB9891", size = 1.2) +
theme_minimal() +
labs(title = "Resultados de la regresión lineal multivariante con selección de modelos",
subtitle = "Los valores predichos deberían estar cercanos a la diagonal",
x = "Predicciones",
y = "Valores reales")
g2 <- pred_test |>
select(.pred, SalePrice) |>
gather(Distribución, value) |>
ggplot(aes(x = value, color = Distribución, fill = Distribución)) +
geom_density(alpha = 0.6) +
theme_minimal() +
labs(title = "Distribución de las predicciones sobre los valores reales de SalePrice",
x = "Distribuciones",
y = "Frecuencia")
multiplot(g1, g2)
Como se puede apreciar en los gráficos, las predicciones son bastante buenas, bastante mejores que en los anteriores modelos. A pesar de que se han incumplido algunos de los supuestos de regresión, la recta se ajusta bastante bien a los datos.
Compararemos ahora los tres modelos de regresión: el modelo univariante, el modelo multivariante saturado y el modelo multivariante con selección de modelos.
compare_performance(reg_fit_uni |> extract_fit_engine(),
reg_fit_multi |> extract_fit_engine(),
reg_fit_bic_test |> extract_fit_engine())
# Comparison of Model Performance Indices
Name | Model | AIC (weights) | AICc (weights) | BIC (weights) | R2 | R2 (adj.) | RMSE | Sigma
------------------------------------------------------------------------------------------------------
..1 | lm | -1932.0 (>.999) | -1931.9 (>.999) | -1916.1 (>.999) | 0.484 | 0.483 | 0.125 | 0.125
..2 | lm | -970.6 (<.001) | -967.4 (<.001) | -798.7 (<.001) | 0.891 | 0.887 | 0.133 | 0.136
..3 | lm | -1030.6 (<.001) | -1023.4 (<.001) | -772.7 (<.001) | 0.902 | 0.896 | 0.126 | 0.130
Como se puede observar, el mejor modelo es la regresión multivariante con selección de modelos. Presenta un \(R^2\) del 0.902 y un \(RMSE\) bastante bajo.
Nos quedamos finalmente con la regresión multivariante con selección de modelos (nos salió un 0.14492 de score :D).
# Predicciones
pred_bic <-
predict(reg_fit_bic, house_test)
summary(pred_bic)
# Visualización de las predicciones para cada `Id`
final_pred <-
data.frame(Id = c(1461:2919), pred_bic) |>
dplyr::rename(SalePrice = .pred) |>
mutate(exp(SalePrice))
head(final_pred)
# Exportamos nuestro dataframe para subirlo a Kaggle
write.csv(final_pred, file = 'house_entrega_bic.csv', row.names = FALSE)